alg_tools: comparison src/linops.rs

-:f75bf34adda0
+:d2994e34a5f5
 use crate::direct_product::Pair;
 use crate::error::DynResult;
 use crate::instance::Instance;
 pub use crate::mapping::{Composition, DifferentiableImpl, Mapping, Space};
-use crate::norms::{Linfinity, Norm, NormExponent, PairNorm, L1, L2};
+use crate::norms::{HasDual, Linfinity, Norm, NormExponent, PairNorm, L1, L2};
 use crate::types::*;
 use numeric_literals::replace_float_literals;
 use serde::Serialize;
 use std::marker::PhantomData;
 use std::ops::Mul;
 //         self.apply(x)
 //     }
 // }
 /// Efficient in-place summation.
-#[replace_float_literals(F::cast_from(literal))]
+#[replace_float_literals(Self::Field::cast_from(literal))]
-pub trait AXPY<F, X = Self>: Space + std::ops::MulAssign<F>
+pub trait AXPY<X = Self>: Space + std::ops::MulAssign<Self::Field>
 where
-F: Num,
 X: Space,
 {
-type Owned: AXPY<F, X>;
+type Field: Num;
+type Owned: AXPY<X, Field = Self::Field>;
+// type OriginGen<'a>: OriginGenerator<Self::Owned, F>
+// where
+//     Self: 'a;
 /// Computes  `y = βy + αx`, where `y` is `Self`.
-fn axpy<I: Instance<X>>(&mut self, α: F, x: I, β: F);
+fn axpy<I: Instance<X>>(&mut self, α: Self::Field, x: I, β: Self::Field);
 /// Copies `x` to `self`.
 fn copy_from<I: Instance<X>>(&mut self, x: I) {
 self.axpy(1.0, x, 0.0)
 }
 /// Computes  `y = αx`, where `y` is `Self`.
-fn scale_from<I: Instance<X>>(&mut self, α: F, x: I) {
+fn scale_from<I: Instance<X>>(&mut self, α: Self::Field, x: I) {
 self.axpy(α, x, 0.0)
 }
 /// Return a similar zero as `self`.
 fn similar_origin(&self) -> Self::Owned;
+// {
+//     self.make_origin_generator().make_origin()
+// }
+/// Return a similar zero as `x`.
+fn similar_origin_inst<I: Instance<Self>>(x: I) -> Self::Owned {
+x.eval(|xr| xr.similar_origin())
+}
 /// Set self to zero.
 fn set_zero(&mut self);
+//fn make_origin_generator(&self) -> Self::OriginGen<'_>;
 }
 /// Efficient in-place application for [`Linear`] operators.
 #[replace_float_literals(F::cast_from(literal))]
 pub trait GEMV<F: Num, X: Space, Y = <Self as Mapping<X>>::Codomain>: Linear<X> {
 impl<X: Clone + Space> Linear<X> for IdOp<X> {}
 #[replace_float_literals(F::cast_from(literal))]
 impl<F: Num, X, Y> GEMV<F, X, Y> for IdOp<X>
 where
-Y: AXPY<F, X>,
+Y: AXPY<X, Field = F>,
 X: Clone + Space,
 {
 // Computes  `y = αAx + βy`, where `A` is `Self`.
 fn gemv<I: Instance<X>>(&self, y: &mut Y, α: F, x: I, β: F) {
 y.axpy(α, x, β)
 fn preadjoint(&self) -> Self::Preadjoint<'_> {
 IdOp::new()
 }
 }
-/// The zero operator
+/// The zero operator from a space to itself
 #[derive(Clone, Copy, Debug, Serialize, Eq, PartialEq)]
-pub struct ZeroOp<'a, X, XD, Y, F> {
+pub struct SimpleZeroOp;
-zero: &'a Y, // TODO: don't pass this in `new`; maybe not even store.
-dual_or_predual_zero: XD,
+impl<X> Mapping<X> for SimpleZeroOp
-_phantoms: PhantomData<(X, Y, F)>,
+where
-}
+X: AXPY + Instance<X>,
+{
-// TODO: Need to make Zero in Instance.
+type Codomain = X::Owned;
-impl<'a, F: Num, X: Space, XD, Y: Space + Clone> ZeroOp<'a, X, XD, Y, F> {
+fn apply<I: Instance<X>>(&self, x: I) -> X::Owned {
-pub fn new(zero: &'a Y, dual_or_predual_zero: XD) -> Self {
+X::similar_origin_inst(x)
-ZeroOp {
+}
-zero,
+}
-dual_or_predual_zero,
-_phantoms: PhantomData,
+impl<X> Linear<X> for SimpleZeroOp where X: AXPY + Instance<X> {}
-}
-}
-}
-impl<'a, F: Num, X: Space, XD, Y: AXPY<F> + Clone> Mapping<X> for ZeroOp<'a, X, XD, Y, F> {
-type Codomain = Y;
-fn apply<I: Instance<X>>(&self, _x: I) -> Y {
-self.zero.clone()
-}
-}
-impl<'a, F: Num, X: Space, XD, Y: AXPY<F> + Clone> Linear<X> for ZeroOp<'a, X, XD, Y, F> {}
 #[replace_float_literals(F::cast_from(literal))]
-impl<'a, F, X, XD, Y> GEMV<F, X, Y> for ZeroOp<'a, X, XD, Y, F>
+impl<X, Y, F> GEMV<F, X, Y> for SimpleZeroOp
 where
 F: Num,
-Y: AXPY<F, Y> + Clone,
+Y: AXPY<Field = F>,
-X: Space,
+X: AXPY<Field = F> + Instance<X>,
 {
 // Computes  `y = αAx + βy`, where `A` is `Self`.
 fn gemv<I: Instance<X>>(&self, y: &mut Y, _α: F, _x: I, β: F) {
 *y *= β;
 }
 fn apply_mut<I: Instance<X>>(&self, y: &mut Y, _x: I) {
 y.set_zero();
 }
 }
-impl<'a, F, X, XD, Y, E1, E2> BoundedLinear<X, E1, E2, F> for ZeroOp<'a, X, XD, Y, F>
+impl<X, F, E1, E2> BoundedLinear<X, E1, E2, F> for SimpleZeroOp
 where
-X: Space + Norm<E1, F>,
+F: Num,
-Y: AXPY<F> + Clone + Norm<E2, F>,
+X: AXPY<Field = F> + Instance<X> + Norm<E1, F>,
-F: Num,
 E1: NormExponent,
 E2: NormExponent,
 {
 fn opnorm_bound(&self, _xexp: E1, _codexp: E2) -> DynResult<F> {
 Ok(F::ZERO)
 }
 }
-impl<'a, F: Num, X, XD, Y, Yprime: Space> Adjointable<X, Yprime> for ZeroOp<'a, X, XD, Y, F>
+impl<X, F> Adjointable<X, X::DualSpace> for SimpleZeroOp
 where
-X: Space,
+F: Num,
-Y: AXPY<F> + Clone + 'static,
+X: AXPY<Field = F> + Instance<X> + HasDual<F>,
-XD: AXPY<F> + Clone + 'static,
+X::DualSpace: AXPY<Owned = X::DualSpace>,
 {
-type AdjointCodomain = XD;
+type AdjointCodomain = X::DualSpace;
 type Adjoint<'b>
-= ZeroOp<'b, Yprime, (), XD, F>
+= SimpleZeroOp
 where
 Self: 'b;
 // () means not (pre)adjointable.
 fn adjoint(&self) -> Self::Adjoint<'_> {
-ZeroOp::new(&self.dual_or_predual_zero, ())
+SimpleZeroOp
 }
 }
-impl<'a, F, X, XD, Y, Ypre> Preadjointable<X, Ypre> for ZeroOp<'a, X, XD, Y, F>
+/*
-where
+/// Trait for forming origins (zeroes) in vector spaces
-F: Num,
+pub trait OriginGenerator<X, F: Num = f64> {
-X: Space,
+fn make_origin(&self) -> X;
-Y: AXPY<F> + Clone,
+}
-Ypre: Space,
-XD: AXPY<F> + Clone + 'static,
+/// Origin generator for statically sized Euclidean spaces.
-{
+pub struct StaticEuclideanOriginGenerator;
-type PreadjointCodomain = XD;
-type Preadjoint<'b>
+impl<X, F: Float> OriginGenerator<X, F> for StaticEuclideanOriginGenerator
-= ZeroOp<'b, Ypre, (), XD, F>
+where
+X: StaticEuclidean<F> + Euclidean<F, Output = X>,
+{
+#[inline]
+fn make_origin(&self) -> X {
+X::origin()
+}
+}
+/// Origin generator arbitrary spaces that implement [`AXPY`], based on a sample vector.
+pub struct SimilarOriginGenerator<'a, X>(&'a X);
+impl<'a, X, F: Float> OriginGenerator<X, F> for SimilarOriginGenerator<'a, X>
+where
+X: AXPY<F, Owned = X>,
+{
+#[inline]
+fn make_origin(&self) -> X {
+self.0.similar_origin()
+}
+}
+pub struct NotAnOriginGenerator;
+/// The zero operator
+#[derive(Clone, Copy, Debug, Serialize, Eq, PartialEq)]
+pub struct ZeroOp<X, Y, YOrigin, XDOrigin = NotAnOriginGenerator, F: Num = f64> {
+y_origin: YOrigin,
+xd_origin: XDOrigin,
+_phantoms: PhantomData<(X, Y, F)>,
+}
+// TODO: Need to make Zero in Instance.
+impl<X, Y, F, YOrigin, XDOrigin> ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+{
+pub fn new(y_origin: YOrigin, xd_origin: XDOrigin) -> Self {
+ZeroOp {
+y_origin,
+xd_origin,
+_phantoms: PhantomData,
+}
+}
+}
+impl<X, Y, F, YOrigin, XDOrigin> Mapping<X> for ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+Y: Space,
+X: Space + Instance<X>,
+{
+type Codomain = Y;
+fn apply<I: Instance<X>>(&self, _: I) -> Y {
+self.y_origin.make_origin()
+}
+}
+impl<X, Y, F, YOrigin, XDOrigin> Linear<X> for ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+Y: Space,
+X: Space + Instance<X>,
+{
+}
+#[replace_float_literals(F::cast_from(literal))]
+impl<X, Y, F, YOrigin, XDOrigin> GEMV<F, X, Y> for ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+Y: Space + AXPY<F>,
+X: Space + Instance<X>,
+{
+// Computes  `y = αAx + βy`, where `A` is `Self`.
+fn gemv<I: Instance<X>>(&self, y: &mut Y, _α: F, _x: I, β: F) {
+*y *= β;
+}
+fn apply_mut<I: Instance<X>>(&self, y: &mut Y, _x: I) {
+y.set_zero();
+}
+}
+impl<X, Y, F, YOrigin, XDOrigin, E1, E2> BoundedLinear<X, E1, E2, F>
+for ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+Y: Space + AXPY<F> + Norm<E2, F>,
+X: Space + Instance<X> + Norm<E1, F>,
+E1: NormExponent,
+E2: NormExponent,
+{
+fn opnorm_bound(&self, _xexp: E1, _codexp: E2) -> DynResult<F> {
+Ok(F::ZERO)
+}
+}
+impl<X, Y, Yprime, F, YOrigin, XDOrigin> Adjointable<X, Yprime>
+for ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+Y: Space + AXPY<F>,
+X: Space + Instance<X> + HasDual<F>,
+XDOrigin: OriginGenerator<X::DualSpace, F> + Clone,
+Yprime: Space + Instance<Yprime>,
+{
+type AdjointCodomain = X::DualSpace;
+type Adjoint<'b>
+= ZeroOp<Yprime, X::DualSpace, XDOrigin, NotAnOriginGenerator, F>
 where
 Self: 'b;
 // () means not (pre)adjointable.
-fn preadjoint(&self) -> Self::Preadjoint<'_> {
+fn adjoint(&self) -> Self::Adjoint<'_> {
-ZeroOp::new(&self.dual_or_predual_zero, ())
+ZeroOp::new(self.xd_origin.clone(), NotAnOriginGenerator)
 }
 }
+*/
+/*
+impl<X, Y, Ypre, Xpre, F, YOrigin, XDOrigin> Preadjointable<X, Ypre>
+for ZeroOp<X, Y, YOrigin, XDOrigin, F>
+where
+F: Num,
+YOrigin: OriginGenerator<Y, F>,
+Y: Space + AXPY<F>,
+X: Space + Instance<X> + HasDual<F>,
+XDOrigin: OriginGenerator<Xpre, F> + Clone,
+Ypre: Space + Instance<Ypre> + HasDual<DualSpace = X>,
+{
+type PreadjointCodomain = Xpre;
+type Preadjoint<'b>
+= ZeroOp<Ypre, Xpre, XDOrigin, NotAnOriginGenerator, F>
+where
+Self: 'b;
+// () means not (pre)adjointable.
+fn adjoint(&self) -> Self::Adjoint<'_> {
+ZeroOp::new(self.xpre_origin.clone(), NotAnOriginGenerator)
+}
+}
+*/
 impl<S, T, E, X> Linear<X> for Composition<S, T, E>
 where
 X: Space,
 T: Linear<X>,

Mercurial > repos > alg_tools / file comparison

comparison: src/linops.rs

src/linops.rs