feanor_math/algorithms/linsolve/
mod.rs

1use std::alloc::{Allocator, Global};
2
3use crate::divisibility::DivisibilityRing;
4use crate::matrix::{AsPointerToSlice, SubmatrixMut};
5use crate::pid::PrincipalIdealRing;
6use crate::ring::*;
7use crate::rings::extension::extension_impl::FreeAlgebraImplBase;
8use crate::seq::VectorView;
9
10use super::convolution::ConvolutionAlgorithm;
11
12///
13/// Contains the algorithm to compute the "pre-smith" form of a matrix, which
14/// is sufficient for solving linear systems. Works over all principal ideal rings.
15/// 
16pub mod smith;
17///
18/// Contains the algorithm for solving linear systems over free ring extensions.
19/// 
20pub mod extension;
21
22///
23/// Result of trying to solve a linear system.
24/// 
25/// Possible values are:
26///  - [`SolveResult::FoundUniqueSolution`]: The system is guaranteed to have a unique solution.
27///  - [`SolveResult::FoundSomeSolution`]: The system has at least one solution. This is also an allowed
28///    value for systems with a unique solution.
29///  - [`SolveResult::NoSolution`]: The system is unsolvable.
30/// 
31#[derive(Debug, Clone, PartialEq, Eq, Hash)]
32pub enum SolveResult {
33    /// The system has at least one solution. This is also an allowed
34    /// value for systems with a unique solution.
35    FoundSomeSolution,
36    /// The system is guaranteed to have a unique solution 
37    FoundUniqueSolution,
38    /// The system has no solution. In particular, the output matrix
39    /// is in an unspecified (but safe) state.
40    NoSolution
41}
42
43impl SolveResult {
44
45    ///
46    /// Returns whether some solution to the system has been found.
47    /// 
48    pub fn is_solved(&self) -> bool {
49        match self {
50            Self::FoundSomeSolution | Self::FoundUniqueSolution => true,
51            Self::NoSolution => false
52        }
53    }
54    
55    ///
56    /// Panics if the system does not have a solution.
57    /// 
58    pub fn assert_solved(&self) {
59        assert!(self.is_solved());
60    }
61}
62
63///
64/// Class for rings over which we can solve linear systems.
65/// 
66pub trait LinSolveRing: DivisibilityRing {
67
68    ///
69    /// Tries to find a matrix `X` such that `lhs * X = rhs`.
70    /// 
71    /// If a solution exists, it will be written to `out`. Otherwise, `out` will have an
72    /// unspecified (but valid) value after the function returns. Similarly, `lhs` and `rhs`
73    /// will be modified in an unspecified way, but its entries will always be valid ring elements.
74    /// 
75    /// Note that if a solution is found, either [`SolveResult::FoundSomeSolution`] or [`SolveResult::FoundUniqueSolution`]
76    /// are returned. If there are multiple solutions, only former will be returned. However, implementations
77    /// are free to also return [`SolveResult::FoundSomeSolution`] in cases where there is a unique solution.
78    /// 
79    fn solve_right<V1, V2, V3, A>(&self, lhs: SubmatrixMut<V1, Self::Element>, rhs: SubmatrixMut<V2, Self::Element>, out: SubmatrixMut<V3, Self::Element>, allocator: A) -> SolveResult
80        where V1: AsPointerToSlice<Self::Element>,
81            V2: AsPointerToSlice<Self::Element>,
82            V3: AsPointerToSlice<Self::Element>,
83            A: Allocator;
84}
85
86///
87/// [`RingStore`] corresponding to [`LinSolveRing`].
88/// 
89pub trait LinSolveRingStore: RingStore
90    where Self::Type: LinSolveRing
91{
92    ///
93    /// Solves a linear system `lhs * X = rhs`.
94    /// 
95    /// For details, see [`LinSolveRing::solve_right()`].
96    /// 
97    fn solve_right<V1, V2, V3>(&self, lhs: SubmatrixMut<V1, El<Self>>, rhs: SubmatrixMut<V2, El<Self>>, out: SubmatrixMut<V3, El<Self>>) -> SolveResult
98        where V1: AsPointerToSlice<El<Self>>,
99            V2: AsPointerToSlice<El<Self>>,
100            V3: AsPointerToSlice<El<Self>>
101    {
102        self.get_ring().solve_right(lhs, rhs, out, Global)
103    }
104
105    ///
106    /// Solves a linear system `lhs * X = rhs`.
107    /// 
108    /// For details, see [`LinSolveRing::solve_right()`].
109    /// 
110    #[stability::unstable(feature = "enable")]
111    fn solve_right_with<V1, V2, V3, A>(&self, lhs: SubmatrixMut<V1, El<Self>>, rhs: SubmatrixMut<V2, El<Self>>, out: SubmatrixMut<V3, El<Self>>, allocator: A) -> SolveResult
112        where V1: AsPointerToSlice<El<Self>>,
113            V2: AsPointerToSlice<El<Self>>,
114            V3: AsPointerToSlice<El<Self>>,
115            A: Allocator
116    {
117        self.get_ring().solve_right(lhs, rhs, out, allocator)
118    }
119}
120
121impl<R> LinSolveRingStore for R
122    where R: RingStore, R::Type: LinSolveRing
123{}
124
125impl<R: ?Sized + PrincipalIdealRing> LinSolveRing for R {
126
127    default fn solve_right<V1, V2, V3, A>(&self, lhs: SubmatrixMut<V1, Self::Element>, rhs: SubmatrixMut<V2, Self::Element>, out: SubmatrixMut<V3, Self::Element>, allocator: A) -> SolveResult
128        where V1: AsPointerToSlice<Self::Element>, V2: AsPointerToSlice<Self::Element>, V3: AsPointerToSlice<Self::Element>,
129            A: Allocator
130    {
131        smith::solve_right_using_pre_smith(RingRef::new(self), lhs, rhs, out, allocator)
132    }
133}
134
135impl<R, V, A_ring, C_ring> LinSolveRing for FreeAlgebraImplBase<R, V, A_ring, C_ring>
136    where R: RingStore,
137        R::Type: LinSolveRing,
138        V: VectorView<El<R>>,
139        A_ring: Allocator + Clone, 
140        C_ring: ConvolutionAlgorithm<R::Type>
141{
142    fn solve_right<V1, V2, V3, A>(&self, lhs: SubmatrixMut<V1, Self::Element>, rhs: SubmatrixMut<V2, Self::Element>, out: SubmatrixMut<V3, Self::Element>, allocator: A) -> SolveResult
143        where V1: AsPointerToSlice<Self::Element>,
144            V2: AsPointerToSlice<Self::Element>,
145            V3: AsPointerToSlice<Self::Element>,
146            A:  Allocator
147    {
148        extension::solve_right_over_extension(RingRef::new(self), lhs, rhs, out, allocator)
149    }
150}
feanor_math/algorithms/linsolve/mod.rs

feanor_math/algorithms/linsolve/
mod.rs
