feanor_math/algorithms/linsolve/
mod.rs

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

feanor_math/algorithms/linsolve/
mod.rs
