Mercurial > repos > alg_tools / file comparison

comparison: src/direct_product.rs

src/direct_product.rs

changeset 90
b3c35d16affe
parent 64
4f6ca107ccb1
child 82
981069ef919b
equal deleted inserted replaced
25:d14c877e14b7 90:b3c35d16affe
1 /*!
2 Direct products of the form $A \times B$.
3
4 TODO: This could be easily much more generic if `derive_more` could derive arithmetic
5 operations on references.
6 */
7
8 use core::ops::{Mul,MulAssign,Div,DivAssign,Add,AddAssign,Sub,SubAssign,Neg};
9 use std::clone::Clone;
10 use serde::{Serialize, Deserialize};
11 use crate::types::{Num, Float};
12 use crate::{maybe_lifetime, maybe_ref};
13 use crate::euclidean::Euclidean;
14 use crate::instance::{Instance, InstanceMut, Decomposition, DecompositionMut, MyCow};
15 use crate::mapping::Space;
16 use crate::linops::AXPY;
17 use crate::loc::Loc;
18 use crate::norms::{Norm, PairNorm, NormExponent, Normed, HasDual, L2};
19
20 #[derive(Debug,Clone,Copy,PartialEq,Eq,Serialize,Deserialize)]
21 pub struct Pair<A, B> (pub A, pub B);
22
23 impl<A, B> Pair<A,B> {
24 pub fn new(a : A, b : B) -> Pair<A,B> { Pair(a, b) }
25 }
26
27 impl<A, B> From<(A,B)> for Pair<A,B> {
28 #[inline]
29 fn from((a, b) : (A, B)) -> Pair<A,B> { Pair(a, b) }
30 }
31
32 impl<A, B> From<Pair<A,B>> for (A,B) {
33 #[inline]
34 fn from(Pair(a, b) : Pair<A, B>) -> (A,B) { (a, b) }
35 }
36
37 macro_rules! impl_binop {
38 (($a : ty, $b : ty), $trait : ident, $fn : ident, $refl:ident, $refr:ident) => {
39 impl_binop!(@doit: $a, $b, $trait, $fn;
40 maybe_lifetime!($refl, &'l Pair<$a,$b>),
41 (maybe_lifetime!($refl, &'l $a),
42 maybe_lifetime!($refl, &'l $b));
43 maybe_lifetime!($refr, &'r Pair<Ai,Bi>),
44 (maybe_lifetime!($refr, &'r Ai),
45 maybe_lifetime!($refr, &'r Bi));
46 $refl, $refr);
47 };
48
49 (@doit: $a:ty, $b:ty,
50 $trait:ident, $fn:ident;
51 $self:ty, ($aself:ty, $bself:ty);
52 $in:ty, ($ain:ty, $bin:ty);
53 $refl:ident, $refr:ident) => {
54 impl<'l, 'r, Ai, Bi> $trait<$in>
55 for $self
56 where $aself: $trait<$ain>,
57 $bself: $trait<$bin> {
58 type Output = Pair<<$aself as $trait<$ain>>::Output,
59 <$bself as $trait<$bin>>::Output>;
60
61 #[inline]
62 fn $fn(self, y : $in) -> Self::Output {
63 Pair(maybe_ref!($refl, self.0).$fn(maybe_ref!($refr, y.0)),
64 maybe_ref!($refl, self.1).$fn(maybe_ref!($refr, y.1)))
65 }
66 }
67 };
68 }
69
70 macro_rules! impl_assignop {
71 (($a : ty, $b : ty), $trait : ident, $fn : ident, $refr:ident) => {
72 impl_assignop!(@doit: $a, $b,
73 $trait, $fn;
74 maybe_lifetime!($refr, &'r Pair<Ai,Bi>),
75 (maybe_lifetime!($refr, &'r Ai),
76 maybe_lifetime!($refr, &'r Bi));
77 $refr);
78 };
79 (@doit: $a : ty, $b : ty,
80 $trait:ident, $fn:ident;
81 $in:ty, ($ain:ty, $bin:ty);
82 $refr:ident) => {
83 impl<'r, Ai, Bi> $trait<$in>
84 for Pair<$a,$b>
85 where $a: $trait<$ain>,
86 $b: $trait<$bin> {
87 #[inline]
88 fn $fn(&mut self, y : $in) -> () {
89 self.0.$fn(maybe_ref!($refr, y.0));
90 self.1.$fn(maybe_ref!($refr, y.1));
91 }
92 }
93 }
94 }
95
96 macro_rules! impl_scalarop {
97 (($a : ty, $b : ty), $field : ty, $trait : ident, $fn : ident, $refl:ident) => {
98 impl_scalarop!(@doit: $field,
99 $trait, $fn;
100 maybe_lifetime!($refl, &'l Pair<$a,$b>),
101 (maybe_lifetime!($refl, &'l $a),
102 maybe_lifetime!($refl, &'l $b));
103 $refl);
104 };
105 (@doit: $field : ty,
106 $trait:ident, $fn:ident;
107 $self:ty, ($aself:ty, $bself:ty);
108 $refl:ident) => {
109 // Scalar as Rhs
110 impl<'l> $trait<$field>
111 for $self
112 where $aself: $trait<$field>,
113 $bself: $trait<$field> {
114 type Output = Pair<<$aself as $trait<$field>>::Output,
115 <$bself as $trait<$field>>::Output>;
116 #[inline]
117 fn $fn(self, a : $field) -> Self::Output {
118 Pair(maybe_ref!($refl, self.0).$fn(a),
119 maybe_ref!($refl, self.1).$fn(a))
120 }
121 }
122 }
123 }
124
125 // Not used due to compiler overflow
126 #[allow(unused_macros)]
127 macro_rules! impl_scalarlhs_op {
128 (($a : ty, $b : ty), $field : ty, $trait:ident, $fn:ident, $refr:ident) => {
129 impl_scalarlhs_op!(@doit: $trait, $fn,
130 maybe_lifetime!($refr, &'r Pair<$a,$b>),
131 (maybe_lifetime!($refr, &'r $a),
132 maybe_lifetime!($refr, &'r $b));
133 $refr, $field);
134 };
135 (@doit: $trait:ident, $fn:ident,
136 $in:ty, ($ain:ty, $bin:ty);
137 $refr:ident, $field:ty) => {
138 impl<'r> $trait<$in>
139 for $field
140 where $field : $trait<$ain>
141 + $trait<$bin> {
142 type Output = Pair<<$field as $trait<$ain>>::Output,
143 <$field as $trait<$bin>>::Output>;
144 #[inline]
145 fn $fn(self, x : $in) -> Self::Output {
146 Pair(self.$fn(maybe_ref!($refr, x.0)),
147 self.$fn(maybe_ref!($refr, x.1)))
148 }
149 }
150 };
151 }
152
153 macro_rules! impl_scalar_assignop {
154 (($a : ty, $b : ty), $field : ty, $trait : ident, $fn : ident) => {
155 impl<'r> $trait<$field>
156 for Pair<$a, $b>
157 where $a: $trait<$field>, $b: $trait<$field> {
158 #[inline]
159 fn $fn(&mut self, a : $field) -> () {
160 self.0.$fn(a);
161 self.1.$fn(a);
162 }
163 }
164 }
165 }
166
167 macro_rules! impl_unaryop {
168 (($a : ty, $b : ty), $trait:ident, $fn:ident, $refl:ident) => {
169 impl_unaryop!(@doit: $trait, $fn;
170 maybe_lifetime!($refl, &'l Pair<$a,$b>),
171 (maybe_lifetime!($refl, &'l $a),
172 maybe_lifetime!($refl, &'l $b));
173 $refl);
174 };
175 (@doit: $trait:ident, $fn:ident;
176 $self:ty, ($aself:ty, $bself:ty);
177 $refl : ident) => {
178 impl<'l> $trait
179 for $self
180 where $aself: $trait,
181 $bself: $trait {
182 type Output = Pair<<$aself as $trait>::Output,
183 <$bself as $trait>::Output>;
184 #[inline]
185 fn $fn(self) -> Self::Output {
186 Pair(maybe_ref!($refl, self.0).$fn(),
187 maybe_ref!($refl, self.1).$fn())
188 }
189 }
190 }
191 }
192
193 #[macro_export]
194 macro_rules! impl_pair_vectorspace_ops {
195 (($a:ty, $b:ty), $field:ty) => {
196 impl_pair_vectorspace_ops!(@binary, ($a, $b), Add, add);
197 impl_pair_vectorspace_ops!(@binary, ($a, $b), Sub, sub);
198 impl_pair_vectorspace_ops!(@assign, ($a, $b), AddAssign, add_assign);
199 impl_pair_vectorspace_ops!(@assign, ($a, $b), SubAssign, sub_assign);
200 impl_pair_vectorspace_ops!(@scalar, ($a, $b), $field, Mul, mul);
201 impl_pair_vectorspace_ops!(@scalar, ($a, $b), $field, Div, div);
202 // Compiler overflow
203 // $(
204 // impl_pair_vectorspace_ops!(@scalar_lhs, ($a, $b), $field, $impl_scalarlhs_op, Mul, mul);
205 // )*
206 impl_pair_vectorspace_ops!(@scalar_assign, ($a, $b), $field, MulAssign, mul_assign);
207 impl_pair_vectorspace_ops!(@scalar_assign, ($a, $b), $field, DivAssign, div_assign);
208 impl_pair_vectorspace_ops!(@unary, ($a, $b), Neg, neg);
209 };
210 (@binary, ($a : ty, $b : ty), $trait : ident, $fn : ident) => {
211 impl_binop!(($a, $b), $trait, $fn, ref, ref);
212 impl_binop!(($a, $b), $trait, $fn, ref, noref);
213 impl_binop!(($a, $b), $trait, $fn, noref, ref);
214 impl_binop!(($a, $b), $trait, $fn, noref, noref);
215 };
216 (@assign, ($a : ty, $b : ty), $trait : ident, $fn :ident) => {
217 impl_assignop!(($a, $b), $trait, $fn, ref);
218 impl_assignop!(($a, $b), $trait, $fn, noref);
219 };
220 (@scalar, ($a : ty, $b : ty), $field : ty, $trait : ident, $fn :ident) => {
221 impl_scalarop!(($a, $b), $field, $trait, $fn, ref);
222 impl_scalarop!(($a, $b), $field, $trait, $fn, noref);
223 };
224 (@scalar_lhs, ($a : ty, $b : ty), $field : ty, $trait : ident, $fn : ident) => {
225 impl_scalarlhs_op!(($a, $b), $field, $trait, $fn, ref);
226 impl_scalarlhs_op!(($a, $b), $field, $trait, $fn, noref);
227 };
228 (@scalar_assign, ($a : ty, $b : ty), $field : ty, $trait : ident, $fn : ident) => {
229 impl_scalar_assignop!(($a, $b), $field, $trait, $fn);
230 };
231 (@unary, ($a : ty, $b : ty), $trait : ident, $fn : ident) => {
232 impl_unaryop!(($a, $b), $trait, $fn, ref);
233 impl_unaryop!(($a, $b), $trait, $fn, noref);
234 };
235 }
236
237 impl_pair_vectorspace_ops!((f32, f32), f32);
238 impl_pair_vectorspace_ops!((f64, f64), f64);
239
240 type PairOutput<F, A, B> = Pair<<A as Euclidean<F>>::Output, <B as Euclidean<F>>::Output>;
241
242 impl<A, B, F> Euclidean<F> for Pair<A, B>
243 where
244 A : Euclidean<F>,
245 B : Euclidean<F>,
246 F : Float,
247 PairOutput<F, A, B> : Euclidean<F>,
248 Self : Sized
249 + Mul<F, Output=PairOutput<F, A, B>> + MulAssign<F>
250 + Div<F, Output=PairOutput<F, A, B>> + DivAssign<F>
251 + Add<Self, Output=PairOutput<F, A, B>>
252 + Sub<Self, Output=PairOutput<F, A, B>>
253 + for<'b> Add<&'b Self, Output=PairOutput<F, A, B>>
254 + for<'b> Sub<&'b Self, Output=PairOutput<F, A, B>>
255 + AddAssign<Self> + for<'b> AddAssign<&'b Self>
256 + SubAssign<Self> + for<'b> SubAssign<&'b Self>
257 + Neg<Output=PairOutput<F, A, B>>
258 {
259 type Output = PairOutput<F, A, B>;
260
261 fn dot<I : Instance<Self>>(&self, other : I) -> F {
262 let Pair(u, v) = other.decompose();
263 self.0.dot(u) + self.1.dot(v)
264 }
265
266 fn norm2_squared(&self) -> F {
267 self.0.norm2_squared() + self.1.norm2_squared()
268 }
269
270 fn dist2_squared<I : Instance<Self>>(&self, other : I) -> F {
271 let Pair(u, v) = other.decompose();
272 self.0.dist2_squared(u) + self.1.dist2_squared(v)
273 }
274 }
275
276 impl<F, A, B, U, V> AXPY<F, Pair<U, V>> for Pair<A, B>
277 where
278 U : Space,
279 V : Space,
280 A : AXPY<F, U>,
281 B : AXPY<F, V>,
282 F : Num,
283 Self : MulAssign<F>,
284 Pair<A, B> : MulAssign<F>,
285 Pair<A::Owned, B::Owned> : AXPY<F, Pair<U, V>>,
286 {
287
288 type Owned = Pair<A::Owned, B::Owned>;
289
290 fn axpy<I : Instance<Pair<U,V>>>(&mut self, α : F, x : I, β : F) {
291 let Pair(u, v) = x.decompose();
292 self.0.axpy(α, u, β);
293 self.1.axpy(α, v, β);
294 }
295
296 fn copy_from<I : Instance<Pair<U,V>>>(&mut self, x : I) {
297 let Pair(u, v) = x.decompose();
298 self.0.copy_from(u);
299 self.1.copy_from(v);
300 }
301
302 fn scale_from<I : Instance<Pair<U,V>>>(&mut self, α : F, x : I) {
303 let Pair(u, v) = x.decompose();
304 self.0.scale_from(α, u);
305 self.1.scale_from(α, v);
306 }
307
308 /// Return a similar zero as `self`.
309 fn similar_origin(&self) -> Self::Owned {
310 Pair(self.0.similar_origin(), self.1.similar_origin())
311 }
312
313 /// Set self to zero.
314 fn set_zero(&mut self) {
315 self.0.set_zero();
316 self.1.set_zero();
317 }
318 }
319
320 /// [`Decomposition`] for working with [`Pair`]s.
321 #[derive(Copy, Clone, Debug)]
322 pub struct PairDecomposition<D, Q>(D, Q);
323
324 impl<A : Space, B : Space> Space for Pair<A, B> {
325 type Decomp = PairDecomposition<A::Decomp, B::Decomp>;
326 }
327
328 impl<A, B, D, Q> Decomposition<Pair<A, B>> for PairDecomposition<D,Q>
329 where
330 A : Space,
331 B : Space,
332 D : Decomposition<A>,
333 Q : Decomposition<B>,
334 {
335 type Decomposition<'b> = Pair<D::Decomposition<'b>, Q::Decomposition<'b>> where Pair<A, B> : 'b;
336 type Reference<'b> = Pair<D::Reference<'b>, Q::Reference<'b>> where Pair<A, B> : 'b;
337
338 #[inline]
339 fn lift<'b>(Pair(u, v) : Self::Reference<'b>) -> Self::Decomposition<'b> {
340 Pair(D::lift(u), Q::lift(v))
341 }
342 }
343
344 impl<A, B, U, V, D, Q> Instance<Pair<A, B>, PairDecomposition<D, Q>> for Pair<U, V>
345 where
346 A : Space,
347 B : Space,
348 D : Decomposition<A>,
349 Q : Decomposition<B>,
350 U : Instance<A, D>,
351 V : Instance<B, Q>,
352 {
353 #[inline]
354 fn decompose<'b>(self)
355 -> <PairDecomposition<D, Q> as Decomposition<Pair<A, B>>>::Decomposition<'b>
356 where Self : 'b, Pair<A, B> : 'b
357 {
358 Pair(self.0.decompose(), self.1.decompose())
359 }
360
361 #[inline]
362 fn ref_instance(&self)
363 -> <PairDecomposition<D, Q> as Decomposition<Pair<A, B>>>::Reference<'_>
364 {
365 Pair(self.0.ref_instance(), self.1.ref_instance())
366 }
367
368 #[inline]
369 fn cow<'b>(self) -> MyCow<'b, Pair<A, B>> where Self : 'b{
370 MyCow::Owned(Pair(self.0.own(), self.1.own()))
371 }
372
373 #[inline]
374 fn own(self) -> Pair<A, B> {
375 Pair(self.0.own(), self.1.own())
376 }
377 }
378
379
380 impl<'a, A, B, U, V, D, Q> Instance<Pair<A, B>, PairDecomposition<D, Q>> for &'a Pair<U, V>
381 where
382 A : Space,
383 B : Space,
384 D : Decomposition<A>,
385 Q : Decomposition<B>,
386 U : Instance<A, D>,
387 V : Instance<B, Q>,
388 &'a U : Instance<A, D>,
389 &'a V : Instance<B, Q>,
390 {
391 #[inline]
392 fn decompose<'b>(self)
393 -> <PairDecomposition<D, Q> as Decomposition<Pair<A, B>>>::Decomposition<'b>
394 where Self : 'b, Pair<A, B> : 'b
395 {
396 Pair(D::lift(self.0.ref_instance()), Q::lift(self.1.ref_instance()))
397 }
398
399 #[inline]
400 fn ref_instance(&self)
401 -> <PairDecomposition<D, Q> as Decomposition<Pair<A, B>>>::Reference<'_>
402 {
403 Pair(self.0.ref_instance(), self.1.ref_instance())
404 }
405
406 #[inline]
407 fn cow<'b>(self) -> MyCow<'b, Pair<A, B>> where Self : 'b {
408 MyCow::Owned(self.own())
409 }
410
411 #[inline]
412 fn own(self) -> Pair<A, B> {
413 let Pair(ref u, ref v) = self;
414 Pair(u.own(), v.own())
415 }
416
417 }
418
419 impl<A, B, D, Q> DecompositionMut<Pair<A, B>> for PairDecomposition<D,Q>
420 where
421 A : Space,
422 B : Space,
423 D : DecompositionMut<A>,
424 Q : DecompositionMut<B>,
425 {
426 type ReferenceMut<'b> = Pair<D::ReferenceMut<'b>, Q::ReferenceMut<'b>> where Pair<A, B> : 'b;
427 }
428
429 impl<A, B, U, V, D, Q> InstanceMut<Pair<A, B>, PairDecomposition<D, Q>> for Pair<U, V>
430 where
431 A : Space,
432 B : Space,
433 D : DecompositionMut<A>,
434 Q : DecompositionMut<B>,
435 U : InstanceMut<A, D>,
436 V : InstanceMut<B, Q>,
437 {
438 #[inline]
439 fn ref_instance_mut(&mut self)
440 -> <PairDecomposition<D, Q> as DecompositionMut<Pair<A, B>>>::ReferenceMut<'_>
441 {
442 Pair(self.0.ref_instance_mut(), self.1.ref_instance_mut())
443 }
444 }
445
446 impl<'a, A, B, U, V, D, Q> InstanceMut<Pair<A, B>, PairDecomposition<D, Q>> for &'a mut Pair<U, V>
447 where
448 A : Space,
449 B : Space,
450 D : DecompositionMut<A>,
451 Q : DecompositionMut<B>,
452 U : InstanceMut<A, D>,
453 V : InstanceMut<B, Q>,
454 {
455 #[inline]
456 fn ref_instance_mut(&mut self)
457 -> <PairDecomposition<D, Q> as DecompositionMut<Pair<A, B>>>::ReferenceMut<'_>
458 {
459 Pair(self.0.ref_instance_mut(), self.1.ref_instance_mut())
460 }
461 }
462
463
464 impl<F, A, B, ExpA, ExpB, ExpJ> Norm<F, PairNorm<ExpA, ExpB, ExpJ>>
465 for Pair<A,B>
466 where
467 F : Num,
468 ExpA : NormExponent,
469 ExpB : NormExponent,
470 ExpJ : NormExponent,
471 A : Norm<F, ExpA>,
472 B : Norm<F, ExpB>,
473 Loc<F, 2> : Norm<F, ExpJ>,
474 {
475 fn norm(&self, PairNorm(expa, expb, expj) : PairNorm<ExpA, ExpB, ExpJ>) -> F {
476 Loc([self.0.norm(expa), self.1.norm(expb)]).norm(expj)
477 }
478 }
479
480
481 impl<F : Float, A, B> Normed<F> for Pair<A,B>
482 where
483 A : Normed<F>,
484 B : Normed<F>,
485 {
486 type NormExp = PairNorm<A::NormExp, B::NormExp, L2>;
487
488 #[inline]
489 fn norm_exponent(&self) -> Self::NormExp {
490 PairNorm(self.0.norm_exponent(), self.1.norm_exponent(), L2)
491 }
492
493 #[inline]
494 fn is_zero(&self) -> bool {
495 self.0.is_zero() && self.1.is_zero()
496 }
497 }
498
499 impl<F : Float, A, B> HasDual<F> for Pair<A,B>
500 where
501 A : HasDual<F>,
502 B : HasDual<F>,
503
504 {
505 type DualSpace = Pair<A::DualSpace, B::DualSpace>;
506 }

Mercurial Repository: alg_tools

mercurial