feanor_math/algorithms/
splitting_field.rs

1use std::alloc::*;
2use std::marker::PhantomData;
3
4use crate::algorithms::convolution::STANDARD_CONVOLUTION;
5use crate::homomorphism::*;
6use crate::matrix::OwnedMatrix;
7use crate::field::*;
8use crate::rings::extension::number_field::*;
9use crate::rings::rational::*;
10use crate::pid::*;
11use crate::primitive_int::StaticRing;
12use crate::rings::extension::extension_impl::*;
13use crate::rings::extension::{FreeAlgebra, FreeAlgebraStore};
14use crate::rings::field::{AsField, AsFieldBase};
15use crate::rings::poly::dense_poly::DensePolyRing;
16use crate::rings::poly::{PolyRing, PolyRingStore};
17use crate::MAX_PROBABILISTIC_REPETITIONS;
18use crate::integer::*;
19use crate::seq::*;
20use crate::ring::*;
21use crate::algorithms::linsolve::*;
22use super::poly_factor::FactorPolyField;
23use super::poly_gcd::PolyTFracGCDRing;
24
25#[stability::unstable(feature = "enable")]
26pub struct NumberFieldHom<R1, Impl1, I1, R2, Impl2, I2>
27    where R1: RingStore<Type = NumberFieldBase<Impl1, I1>>,
28        Impl1: RingStore,
29        Impl1::Type: Field + FreeAlgebra,
30        <Impl1::Type as RingExtension>::BaseRing: RingStore<Type = RationalFieldBase<I1>>,
31        I1: RingStore,
32        I1::Type: IntegerRing,
33        R2: RingStore<Type = NumberFieldBase<Impl2, I2>>,
34        Impl2: RingStore,
35        Impl2::Type: Field + FreeAlgebra,
36        <Impl2::Type as RingExtension>::BaseRing: RingStore<Type = RationalFieldBase<I2>>,
37        I2: RingStore,
38        I2::Type: IntegerRing
39{
40    from: R1,
41    to: R2,
42    from_int: PhantomData<I1>,
43    to_int: PhantomData<I2>,
44    generator_image: El<NumberField<Impl2, I2>>,
45}
46
47impl<R1, Impl1, I1, R2, Impl2, I2> Homomorphism<NumberFieldBase<Impl1, I1>, NumberFieldBase<Impl2, I2>> for NumberFieldHom<R1, Impl1, I1, R2, Impl2, I2>
48    where R1: RingStore<Type = NumberFieldBase<Impl1, I1>>,
49        Impl1: RingStore,
50        Impl1::Type: Field + FreeAlgebra,
51        <Impl1::Type as RingExtension>::BaseRing: RingStore<Type = RationalFieldBase<I1>>,
52        I1: RingStore,
53        I1::Type: IntegerRing,
54        R2: RingStore<Type = NumberFieldBase<Impl2, I2>>,
55        Impl2: RingStore,
56        Impl2::Type: Field + FreeAlgebra,
57        <Impl2::Type as RingExtension>::BaseRing: RingStore<Type = RationalFieldBase<I2>>,
58        I2: RingStore,
59        I2::Type: IntegerRing
60{
61    type DomainStore = R1;
62    type CodomainStore = R2;
63
64    fn domain<'a>(&'a self) -> &'a Self::DomainStore {
65        &self.from
66    }
67
68    fn codomain<'a>(&'a self) -> &'a Self::CodomainStore {
69        &self.to
70    }
71
72    fn map(&self, x: <NumberFieldBase<Impl1, I1> as RingBase>::Element) -> <NumberFieldBase<Impl2, I2> as RingBase>::Element {
73        self.map_ref(&x)
74    }
75
76    fn map_ref(&self, x: &<NumberFieldBase<Impl1, I1> as RingBase>::Element) -> <NumberFieldBase<Impl2, I2> as RingBase>::Element {
77        let poly_ring = DensePolyRing::new(self.to.base_ring(), "X");
78        return poly_ring.evaluate(
79            &self.from.poly_repr(&poly_ring, &x, self.to.base_ring().can_hom(self.from.base_ring()).unwrap()),
80            &self.generator_image,
81            self.to.inclusion()
82        )
83    }
84}
85
86///
87/// Constructs the field `F[X]/(f(X))` that is isomorphic to `(F[X]/(g(X))[Y]/(h(Y))`
88/// where `F[X]/(g(X))` is the base ring of `poly_ring` and `h` is the given irreducible
89/// polynomial over `F[X]/(g(X))`. 
90/// 
91#[stability::unstable(feature = "enable")]
92pub fn extend_number_field<P>(poly_ring: P, irred_poly: &El<P>) -> (
93    NumberFieldHom<<<P as RingStore>::Type as RingExtension>::BaseRing, DefaultNumberFieldImpl, BigIntRing, NumberField, DefaultNumberFieldImpl, BigIntRing>,
94    El<NumberField>
95)
96    where P: RingStore,
97        P::Type: PolyRing + EuclideanRing,
98        <P::Type as RingExtension>::BaseRing: Clone + RingStore<Type = NumberFieldBase<DefaultNumberFieldImpl, BigIntRing>>
99{
100    assert!(!poly_ring.is_zero(&irred_poly));
101    assert!(poly_ring.degree(&irred_poly).unwrap() > 1);
102    assert!(<NumberFieldBase<_, _> as FactorPolyField>::is_irred(&poly_ring, irred_poly));
103
104    extend_number_field_promise_is_irreducible(poly_ring, irred_poly)
105}
106
107///
108/// Given a number field `K` and an irreducible polynomial `f`, computes a representation of
109/// the number field `L = K[X]/(f)`. The result is returned by the inclusion `K -> L` and 
110/// the element that corresponds to the coset of `X`, i.e. a root of `f` in `L`. Note that the 
111/// canonical generator of `L` does not have to be a root of `f` (this might even be impossible,
112/// e.g. if `f in ZZ[X]` but `K != QQ`).
113/// 
114#[stability::unstable(feature = "enable")]
115pub fn extend_number_field_promise_is_irreducible<P>(poly_ring: P, irred_poly: &El<P>) -> (
116    NumberFieldHom<<<P as RingStore>::Type as RingExtension>::BaseRing, DefaultNumberFieldImpl, BigIntRing, NumberField, DefaultNumberFieldImpl, BigIntRing>,
117    El<NumberField>
118)
119    where P: RingStore,
120        P::Type: PolyRing + EuclideanRing,
121        <P::Type as RingExtension>::BaseRing: Clone + RingStore<Type = NumberFieldBase<DefaultNumberFieldImpl, BigIntRing>>
122{
123    static_assert_impls!(NumberFieldBase<DefaultNumberFieldImpl, BigIntRing>: FactorPolyField);
124    static_assert_impls!(NumberFieldBase<DefaultNumberFieldImpl, BigIntRing>: PolyTFracGCDRing);
125
126    assert!(!poly_ring.is_zero(&irred_poly));
127    assert!(poly_ring.degree(&irred_poly).unwrap() > 1);
128
129    let K = poly_ring.base_ring();
130    let QQ = K.base_ring();
131
132    let L = AsField::from(
133        AsFieldBase::promise_is_perfect_field(
134            FreeAlgebraImpl::new_with(
135                K, 
136                poly_ring.degree(&irred_poly).unwrap(), 
137                (0..poly_ring.degree(&irred_poly).unwrap()).map(|i| K.negate(K.clone_el(poly_ring.coefficient_at(&irred_poly, i)))).collect::<Vec<_>>(),
138                "X",
139                Global,
140                STANDARD_CONVOLUTION
141            )
142        )
143    );
144    
145    let total_rank = K.rank() * L.rank();
146
147    // the main task is to find a primitive element of `extension_ring`, i.e. that generates it over `base_ring`.
148    // we use the following approach:
149    //  - Consider generator `a` of `K` and `b` of `L` over `K`
150    //  - Consider the set `A_n = { b, b + a, b + 2a, ..., b + (n - 1)a }`
151    //  - Consider also the subset `B_n = { x in A_n | QQ[x] != L }`
152    //  - The idea is now to choose a random element from `A_n`, and hope that it is not in `B_n`
153    //  - To estimate the probability, consider the map `B_n -> { maximal proper subfields of L }` that maps any `x`
154    //    to some maximal proper subfield containing `QQ[x]`
155    //  - Then this map is injective, as any field containing `a + ib` and `a + jb`, `j > i` must contain `a, b`, thus be `L`
156    //  - In other words, to find `n`, we need a bound on the maximal proper subfields
157    //  - I believe the degree `[L : QQ]` is such a bound (in the Galois case it is, at least)
158
159    // take `A` twice as large, so that we find a good element with probability >= 1/2
160    let size_of_A = (2 * total_rank) as i32;
161
162    let mut rng = oorandom::Rand64::new(1);
163
164    for _ in 0..MAX_PROBABILISTIC_REPETITIONS {
165
166        let a = StaticRing::<i32>::RING.get_uniformly_random(&size_of_A, || rng.rand_u64());
167        let potential_primitive_element = L.add(L.canonical_gen(), L.inclusion().map(K.int_hom().mul_map(K.canonical_gen(), a)));
168    
169        // lhs is the matrix whose columns are the powers of `b + ka`, w.r.t. the tensor basis
170        let mut lhs = OwnedMatrix::zero(total_rank, total_rank, QQ);
171        let mut current = L.one();
172        for j in 0..total_rank {
173            let current_wrt_basis = L.wrt_canonical_basis(&current);
174            for i1 in 0..L.rank() {
175                let c = current_wrt_basis.at(i1);
176                let c_wrt_basis = K.wrt_canonical_basis(&c);
177                for i2 in 0..K.rank() {
178                    *lhs.at_mut(i1 * K.rank() + i2, j) = c_wrt_basis.at(i2);
179                }
180            }
181            L.mul_assign_ref(&mut current, &potential_primitive_element);
182        }
183    
184        // rhs has three columns: 
185        // - first the `total_rank`-th power of `x + ka` -> this will later give us the minpoly
186        // - second the generator of `K` -> this will give us a repr of this generator in the result field
187        // - third the generator of `L` -> this will give us a repr of this generator in the result field
188        let mut rhs = OwnedMatrix::zero(total_rank, 3, QQ);
189        let current_wrt_basis = L.wrt_canonical_basis(&current);
190        for i1 in 0..L.rank() {
191            let c = current_wrt_basis.at(i1);
192            let c_wrt_basis = K.wrt_canonical_basis(&c);
193            for i2 in 0..K.rank() {
194                *rhs.at_mut(i1 * K.rank() + i2, 0) = c_wrt_basis.at(i2);
195            }
196        }
197        if K.rank() > 1 {
198            *rhs.at_mut(1, 1) = QQ.one();
199        } else {
200            *rhs.at_mut(0, 1) = K.wrt_canonical_basis(&K.canonical_gen()).at(0);
201        }
202        *rhs.at_mut(K.rank(), 2) = QQ.one();
203    
204        let mut sol = OwnedMatrix::zero(total_rank, 3, QQ);
205        let has_sol = QQ.solve_right(lhs.data_mut(), rhs.data_mut(), sol.data_mut()).is_solved();
206
207        if has_sol {
208    
209            let QQX = DensePolyRing::new(QQ, "X");
210            let gen_poly = QQX.from_terms((0..total_rank).map(|i| (QQ.negate(QQ.clone_el(sol.at(i, 0))), i)).chain([(QQ.one(), total_rank)].into_iter()));
211    
212            if let Some((result, x)) = NumberField::try_adjoin_root(&QQX, &gen_poly) {
213
214                // note that `sol` contains coefficients w.r.t. `x` and not `result.canonical_gen()`
215                let K_generator = <_ as RingStore>::sum(&result, (0..total_rank).map(|i| result.inclusion().mul_map(result.pow(result.clone_el(&x), i), QQ.clone_el(sol.at(i, 1)))));
216                let L_generator = <_ as RingStore>::sum(&result, (0..total_rank).map(|i| result.inclusion().mul_map(result.pow(result.clone_el(&x), i), QQ.clone_el(sol.at(i, 2)))));
217
218                let result = NumberFieldHom {
219                    from: K.clone(),
220                    to: result,
221                    from_int: PhantomData,
222                    to_int: PhantomData,
223                    generator_image: K_generator
224                };
225
226                return (result, L_generator);
227            }
228        }
229    }
230
231    unreachable!()
232}
233
234#[cfg(test)]
235use crate::wrapper::RingElementWrapper;
236
237#[test]
238fn test_extend_field() {
239    let ZZ = BigIntRing::RING;
240    let QQ = RationalField::new(ZZ);
241    let ZZX = DensePolyRing::new(&ZZ, "X");
242    let [f] = ZZX.with_wrapped_indeterminate(|X| [X.pow_ref(2) + 1]);
243    let K = NumberField::new(ZZX, &f);
244    let KX = DensePolyRing::new(&K, "X");
245
246    // extend `QQ[i]` by `X^4 - i`
247    let [g] = KX.with_wrapped_indeterminate(|X| [X.pow_ref(4) - RingElementWrapper::new(&KX, KX.inclusion().map(K.canonical_gen()))]);
248    let (extension_field_embedding, x) = extend_number_field(&KX, &g);
249    let ext_field = extension_field_embedding.codomain();
250    assert_eq!(8, ext_field.rank());
251    assert_el_eq!(ext_field, ext_field.neg_one(), ext_field.pow(extension_field_embedding.map(K.canonical_gen()), 2));
252    assert_el_eq!(ext_field, extension_field_embedding.map(K.canonical_gen()), ext_field.pow(x, 4));
253
254    crate::homomorphism::generic_tests::test_homomorphism_axioms(
255        &extension_field_embedding, 
256        [(0, 0), (0, 1), (1, 0), (2, 0), (-1, 0), (0, -1), (-1, -1)].into_iter().map(|(a, b)| K.from_canonical_basis([QQ.int_hom().map(a), QQ.int_hom().map(b)]))
257    );
258
259    // extend `QQ[i]` by `X^3 - 2`
260    let [g] = KX.with_wrapped_indeterminate(|X| [X.pow_ref(3) - 2]);
261    let (extension_field_embedding, x) = extend_number_field(&KX, &g);
262    let ext_field = extension_field_embedding.codomain();
263    assert_eq!(6, ext_field.rank());
264    assert_el_eq!(ext_field, ext_field.neg_one(), ext_field.pow(extension_field_embedding.map(K.canonical_gen()), 2));
265    assert_el_eq!(ext_field, ext_field.int_hom().map(2), ext_field.pow(x, 3));
266    
267    crate::homomorphism::generic_tests::test_homomorphism_axioms(
268        &extension_field_embedding, 
269        [(0, 0), (0, 1), (1, 0), (2, 0), (-1, 0), (0, -1), (-1, -1)].into_iter().map(|(a, b)| K.from_canonical_basis([QQ.int_hom().map(a), QQ.int_hom().map(b)]))
270    );
271}