rlst/operator/
abstract_operator.rs

1//! Definition of a general linear operator.
2
3use crate::dense::types::RlstResult;
4use crate::operator::{FieldType, LinearSpace};
5use num::{One, Zero};
6use std::fmt::Debug;
7
8/// A base operator trait.
9pub trait OperatorBase: Debug {
10    /// Domain space type
11    type Domain: LinearSpace;
12    /// Range space type
13    type Range: LinearSpace;
14
15    /// Get the domain
16    fn domain(&self) -> &Self::Domain;
17    /// Get the range
18    fn range(&self) -> &Self::Range;
19    /// Convert to RLST reference
20    fn as_ref_obj(&self) -> RlstOperatorReference<'_, Self>
21    where
22        Self: Sized,
23    {
24        RlstOperatorReference(self)
25    }
26    /// Form a new operator alpha * self.
27    fn scale(self, alpha: <Self::Range as LinearSpace>::F) -> ScalarTimesOperator<Self>
28    where
29        Self: Sized,
30    {
31        ScalarTimesOperator(self, alpha)
32    }
33    /// Form a new operator self + other.
34    fn sum<Op: OperatorBase<Domain = Self::Domain, Range = Self::Range> + Sized>(
35        self,
36        other: Op,
37    ) -> OperatorSum<Self::Domain, Self::Range, Self, Op>
38    where
39        Self: Sized,
40    {
41        OperatorSum(self, other)
42    }
43
44    /// Form a new operator self * other.
45    fn product<Op: OperatorBase<Range = Self::Domain>>(
46        self,
47        other: Op,
48    ) -> OperatorProduct<Self::Domain, Op, Self>
49    where
50        Self: Sized,
51    {
52        OperatorProduct {
53            op1: other,
54            op2: self,
55        }
56    }
57}
58
59/// Apply an operator as y -> alpha * Ax + beta y
60pub trait AsApply: OperatorBase {
61    /// Apply an operator as y -> alpha * Ax + beta y
62    fn apply_extended(
63        &self,
64        alpha: <Self::Range as LinearSpace>::F,
65        x: &<Self::Domain as LinearSpace>::E,
66        beta: <Self::Range as LinearSpace>::F,
67        y: &mut <Self::Range as LinearSpace>::E,
68    ) -> RlstResult<()>;
69
70    /// Apply
71    fn apply(&self, x: &<Self::Domain as LinearSpace>::E) -> <Self::Range as LinearSpace>::E {
72        let mut out = self.range().zero();
73
74        self.apply_extended(
75            <FieldType<Self::Range> as One>::one(),
76            x,
77            <FieldType<Self::Range> as Zero>::zero(),
78            &mut out,
79        )
80        .unwrap();
81        out
82    }
83}
84
85/// Operator reference
86pub struct RlstOperatorReference<'a, Op: OperatorBase>(&'a Op);
87
88impl<Op: OperatorBase> std::fmt::Debug for RlstOperatorReference<'_, Op> {
89    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
90        f.debug_tuple("OperatorReference").field(&self.0).finish()
91    }
92}
93
94impl<Op: OperatorBase> OperatorBase for RlstOperatorReference<'_, Op> {
95    type Domain = Op::Domain;
96
97    type Range = Op::Range;
98
99    fn domain(&self) -> &Self::Domain {
100        self.0.domain()
101    }
102
103    fn range(&self) -> &Self::Range {
104        self.0.range()
105    }
106}
107
108impl<Op: AsApply> AsApply for RlstOperatorReference<'_, Op> {
109    fn apply_extended(
110        &self,
111        alpha: <Self::Range as LinearSpace>::F,
112        x: &<Self::Domain as LinearSpace>::E,
113        beta: <Self::Range as LinearSpace>::F,
114        y: &mut <Self::Range as LinearSpace>::E,
115    ) -> RlstResult<()> {
116        self.0.apply_extended(alpha, x, beta, y)
117    }
118}
119
120/// Operator sum
121pub struct OperatorSum<
122    Domain: LinearSpace,
123    Range: LinearSpace,
124    Op1: OperatorBase<Domain = Domain, Range = Range>,
125    Op2: OperatorBase<Domain = Domain, Range = Range>,
126>(Op1, Op2);
127
128impl<
129        Domain: LinearSpace,
130        Range: LinearSpace,
131        Op1: OperatorBase<Domain = Domain, Range = Range>,
132        Op2: OperatorBase<Domain = Domain, Range = Range>,
133    > std::fmt::Debug for OperatorSum<Domain, Range, Op1, Op2>
134{
135    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
136        f.debug_tuple("OperatorSum")
137            .field(&&self.0)
138            .field(&&self.1)
139            .finish()
140    }
141}
142
143impl<
144        Domain: LinearSpace,
145        Range: LinearSpace,
146        Op1: OperatorBase<Domain = Domain, Range = Range>,
147        Op2: OperatorBase<Domain = Domain, Range = Range>,
148    > OperatorBase for OperatorSum<Domain, Range, Op1, Op2>
149{
150    type Domain = Domain;
151
152    type Range = Range;
153
154    fn domain(&self) -> &Self::Domain {
155        self.0.domain()
156    }
157
158    fn range(&self) -> &Self::Range {
159        self.0.range()
160    }
161}
162
163impl<
164        Domain: LinearSpace,
165        Range: LinearSpace,
166        Op1: AsApply<Domain = Domain, Range = Range>,
167        Op2: AsApply<Domain = Domain, Range = Range>,
168    > AsApply for OperatorSum<Domain, Range, Op1, Op2>
169{
170    fn apply_extended(
171        &self,
172        alpha: <Self::Range as LinearSpace>::F,
173        x: &<Self::Domain as LinearSpace>::E,
174        beta: <Self::Range as LinearSpace>::F,
175        y: &mut <Self::Range as LinearSpace>::E,
176    ) -> RlstResult<()> {
177        self.0.apply_extended(alpha, x, beta, y)?;
178        self.1
179            .apply_extended(alpha, x, <<Self::Range as LinearSpace>::F as One>::one(), y)?;
180        Ok(())
181    }
182}
183
184/// Operator muiltiplied by a scalar
185pub struct ScalarTimesOperator<Op: OperatorBase>(Op, <Op::Range as LinearSpace>::F);
186
187impl<Op: OperatorBase> std::fmt::Debug for ScalarTimesOperator<Op> {
188    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
189        f.debug_tuple("OperatorSum")
190            .field(&self.0)
191            .field(&self.1)
192            .finish()
193    }
194}
195
196impl<Op: OperatorBase> OperatorBase for ScalarTimesOperator<Op> {
197    type Domain = Op::Domain;
198
199    type Range = Op::Range;
200
201    fn domain(&self) -> &Self::Domain {
202        self.0.domain()
203    }
204
205    fn range(&self) -> &Self::Range {
206        self.0.range()
207    }
208}
209
210impl<Op: AsApply> AsApply for ScalarTimesOperator<Op> {
211    fn apply_extended(
212        &self,
213        alpha: <Self::Range as LinearSpace>::F,
214        x: &<Self::Domain as LinearSpace>::E,
215        beta: <Self::Range as LinearSpace>::F,
216        y: &mut <Self::Range as LinearSpace>::E,
217    ) -> RlstResult<()> {
218        self.0.apply_extended(self.1 * alpha, x, beta, y)?;
219        Ok(())
220    }
221}
222
223/// The product op2 * op1.
224pub struct OperatorProduct<
225    Space: LinearSpace,
226    Op1: OperatorBase<Range = Space>,
227    Op2: OperatorBase<Domain = Space>,
228> {
229    op1: Op1,
230    op2: Op2,
231}
232
233impl<Space: LinearSpace, Op1: OperatorBase<Range = Space>, Op2: OperatorBase<Domain = Space>> Debug
234    for OperatorProduct<Space, Op1, Op2>
235{
236    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
237        f.debug_struct("ProductOperator")
238            .field("op1", &self.op1)
239            .field("op2", &self.op2)
240            .finish()
241    }
242}
243
244impl<Space: LinearSpace, Op1: OperatorBase<Range = Space>, Op2: OperatorBase<Domain = Space>>
245    OperatorBase for OperatorProduct<Space, Op1, Op2>
246{
247    type Domain = Op1::Domain;
248
249    type Range = Op2::Range;
250
251    fn domain(&self) -> &Self::Domain {
252        self.op1.domain()
253    }
254
255    fn range(&self) -> &Self::Range {
256        self.op2.range()
257    }
258}
259
260impl<Space: LinearSpace, Op1: AsApply<Range = Space>, Op2: AsApply<Domain = Space>> AsApply
261    for OperatorProduct<Space, Op1, Op2>
262{
263    fn apply_extended(
264        &self,
265        alpha: <Self::Range as LinearSpace>::F,
266        x: &<Self::Domain as LinearSpace>::E,
267        beta: <Self::Range as LinearSpace>::F,
268        y: &mut <Self::Range as LinearSpace>::E,
269    ) -> RlstResult<()> {
270        self.op2.apply_extended(alpha, &self.op1.apply(x), beta, y)
271    }
272}