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    FoundSomeSolution, FoundUniqueSolution, NoSolution
34}
35
36impl SolveResult {
37
38    pub fn is_solved(&self) -> bool {
39        match self {
40            Self::FoundSomeSolution | Self::FoundUniqueSolution => true,
41            Self::NoSolution => false
42        }
43    }
44    
45    pub fn assert_solved(&self) {
46        assert!(self.is_solved());
47    }
48}
49
50///
51/// Class for rings over which we can solve linear systems.
52/// 
53pub trait LinSolveRing: DivisibilityRing {
54
55    ///
56    /// Tries to find a matrix `X` such that `lhs * X = rhs`.
57    /// 
58    /// If a solution exists, it will be written to `out`. Otherwise, `out` will have an
59    /// unspecified (but valid) value after the function returns. Similarly, `lhs` and `rhs`
60    /// will be modified in an unspecified way, but its entries will always be valid ring elements.
61    /// 
62    /// Note that if a solution is found, either [`SolveResult::FoundSomeSolution`] or [`SolveResult::FoundUniqueSolution`]
63    /// are returned. If there are multiple solutions, only former will be returned. However, implementations
64    /// are free to also return [`SolveResult::FoundSomeSolution`] in cases where there is a unique solution.
65    /// 
66    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
67        where V1: AsPointerToSlice<Self::Element>,
68            V2: AsPointerToSlice<Self::Element>,
69            V3: AsPointerToSlice<Self::Element>,
70            A: Allocator;
71}
72
73///
74/// [`RingStore`] corresponding to [`LinSolveRing`].
75/// 
76pub trait LinSolveRingStore: RingStore
77    where Self::Type: LinSolveRing
78{
79    ///
80    /// Solves a linear system `lhs * X = rhs`.
81    /// 
82    /// For details, see [`LinSolveRing::solve_right()`].
83    /// 
84    fn solve_right<V1, V2, V3>(&self, lhs: SubmatrixMut<V1, El<Self>>, rhs: SubmatrixMut<V2, El<Self>>, out: SubmatrixMut<V3, El<Self>>) -> SolveResult
85        where V1: AsPointerToSlice<El<Self>>,
86            V2: AsPointerToSlice<El<Self>>,
87            V3: AsPointerToSlice<El<Self>>
88    {
89        self.get_ring().solve_right(lhs, rhs, out, Global)
90    }
91
92    ///
93    /// Solves a linear system `lhs * X = rhs`.
94    /// 
95    /// For details, see [`LinSolveRing::solve_right()`].
96    /// 
97    #[stability::unstable(feature = "enable")]
98    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
99        where V1: AsPointerToSlice<El<Self>>,
100            V2: AsPointerToSlice<El<Self>>,
101            V3: AsPointerToSlice<El<Self>>,
102            A: Allocator
103    {
104        self.get_ring().solve_right(lhs, rhs, out, allocator)
105    }
106}
107
108impl<R> LinSolveRingStore for R
109    where R: RingStore, R::Type: LinSolveRing
110{}
111
112impl<R: ?Sized + PrincipalIdealRing> LinSolveRing for R {
113
114    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
115        where V1: AsPointerToSlice<Self::Element>, V2: AsPointerToSlice<Self::Element>, V3: AsPointerToSlice<Self::Element>,
116            A: Allocator
117    {
118        smith::solve_right_using_pre_smith(RingRef::new(self), lhs, rhs, out, allocator)
119    }
120}
121
122impl<R, V, A_ring, C_ring> LinSolveRing for FreeAlgebraImplBase<R, V, A_ring, C_ring>
123    where R: RingStore,
124        R::Type: LinSolveRing,
125        V: VectorView<El<R>>,
126        A_ring: Allocator + Clone, 
127        C_ring: ConvolutionAlgorithm<R::Type>
128{
129    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
130        where V1: AsPointerToSlice<Self::Element>,
131            V2: AsPointerToSlice<Self::Element>,
132            V3: AsPointerToSlice<Self::Element>,
133            A:  Allocator
134    {
135        extension::solve_right_over_extension(RingRef::new(self), lhs, rhs, out, allocator)
136    }
137}
feanor_math/algorithms/linsolve/mod.rs

feanor_math/algorithms/linsolve/
mod.rs
