Givaro Field Elements
Sage includes the Givaro finite field library, for highly optimized arithmetic in finite fields.
Note
The arithmetic is performed by the Givaro C++ library which uses Zech
logs internally to represent finite field elements. This
implementation is the default finite extension field implementation in
Sage for the cardinality less than , as it is vastly faster than
the PARI implementation which uses polynomials to represent finite field
elements. Some functionality in this class however is implemented
using the PARI implementation.
EXAMPLES:
sage: k = GF(5); type(k)
<class 'sage.rings.finite_rings.finite_field_prime_modn.FiniteField_prime_modn_with_category'>
sage: k = GF(5^2,'c'); type(k)
<class 'sage.rings.finite_rings.finite_field_givaro.FiniteField_givaro_with_category'>
sage: k = GF(2^16,'c'); type(k)
<class 'sage.rings.finite_rings.finite_field_ntl_gf2e.FiniteField_ntl_gf2e_with_category'>
sage: k = GF(3^16,'c'); type(k)
<class 'sage.rings.finite_rings.finite_field_pari_ffelt.FiniteField_pari_ffelt_with_category'>
sage: n = previous_prime_power(2^16 - 1)
sage: while is_prime(n):
... n = previous_prime_power(n)
sage: factor(n)
251^2
sage: k = GF(n,'c'); type(k)
<class 'sage.rings.finite_rings.finite_field_givaro.FiniteField_givaro_with_category'>
AUTHORS:
Bases: sage.structure.sage_object.SageObject
Finite Field.
These are implemented using Zech logs and the
cardinality must be less than . By default conway polynomials
are used as minimal polynomial.
INPUT:
q – (must be prime power)
name – variable used for poly_repr (default: 'a')
modulus – you may provide a polynomial to use for reduction or one of the following strings:
Furthermore, for binary fields we allow two more options:
repr – (default: ‘poly’) controls the way elements are printed to the user:
cache – (default: False) if True a cache of all elements of this field is created. Thus, arithmetic does not create new elements which speeds calculations up. Also, if many elements are needed during a calculation this cache reduces the memory requirement as at most order() elements are created.
OUTPUT:
Givaro finite field with characteristic and cardinality
.
EXAMPLES:
By default Conway polynomials are used:
sage: k.<a> = GF(2**8)
sage: -a ^ k.degree()
a^4 + a^3 + a^2 + 1
sage: f = k.modulus(); f
x^8 + x^4 + x^3 + x^2 + 1
You may enforce a modulus:
sage: P.<x> = PolynomialRing(GF(2))
sage: f = x^8 + x^4 + x^3 + x + 1 # Rijndael polynomial
sage: k.<a> = GF(2^8, modulus=f)
sage: k.modulus()
x^8 + x^4 + x^3 + x + 1
sage: a^(2^8)
a
You may enforce a random modulus:
sage: k = GF(3**5, 'a', modulus='random')
sage: k.modulus() # random polynomial
x^5 + 2*x^4 + 2*x^3 + x^2 + 2
For binary fields, you may ask for a minimal weight polynomial:
sage: k = GF(2**10, 'a', modulus='minimal_weight')
sage: k.modulus()
x^10 + x^3 + 1
Three different representations are possible:
sage: sage.rings.finite_rings.finite_field_givaro.FiniteField_givaro(9,repr='poly').gen()
a
sage: sage.rings.finite_rings.finite_field_givaro.FiniteField_givaro(9,repr='int').gen()
3
sage: sage.rings.finite_rings.finite_field_givaro.FiniteField_givaro(9,repr='log').gen()
1
Return a*b - c.
INPUT:
EXAMPLES:
sage: k.<a> = GF(3**3)
sage: k._cache.a_times_b_minus_c(a,a,k(1))
a^2 + 2
Return a*b + c. This is faster than multiplying a and b first and adding c to the result.
INPUT:
EXAMPLES:
sage: k.<a> = GF(2**8)
sage: k._cache.a_times_b_plus_c(a,a,k(1))
a^2 + 1
Return c - a*b.
INPUT:
EXAMPLES:
sage: k.<a> = GF(3**3)
sage: k._cache.c_minus_a_times_b(a,a,k(1))
2*a^2 + 1
Return the characteristic of this field.
EXAMPLES:
sage: p = GF(19^3,'a')._cache.characteristic(); p
19
Coerces several data types to self.
INPUT:
EXAMPLES:
sage: k = GF(2**8, 'a')
sage: e = k.vector_space().gen(1); e
(0, 1, 0, 0, 0, 0, 0, 0)
sage: k(e) #indirect doctest
a
For more examples, see finite_field_givaro.FiniteField_givaro._element_constructor_
Returns the degree of this field over .
EXAMPLES:
sage: K.<a> = GF(9); K._cache.exponent()
2
Given an integer n return a finite field element in self which equals n under the condition that gen() is set to characteristic().
EXAMPLES:
sage: k.<a> = GF(2^8)
sage: k._cache.fetch_int(8)
a^3
sage: e = k._cache.fetch_int(151); e
a^7 + a^4 + a^2 + a + 1
sage: 2^7 + 2^4 + 2^2 + 2 + 1
151
Returns a generator of the field.
EXAMPLES:
sage: K.<a> = GF(625)
sage: K._cache.gen()
a
Given an integer this method returns
where
satisfies
where
is the generator and
is the
characteristic of self.
INPUT:
OUTPUT:
log representation of n
EXAMPLES:
sage: k = GF(7**3, 'a')
sage: k._cache.int_to_log(4)
228
sage: k._cache.int_to_log(3)
57
sage: k.gen()^57
3
Given an integer this method returns
where
satisfies
where
is the generator of self; the
result is interpreted as an integer.
INPUT:
OUTPUT:
integer representation of a finite field element.
EXAMPLES:
sage: k = GF(2**8, 'a')
sage: k._cache.log_to_int(4)
16
sage: k._cache.log_to_int(20)
180
Returns the order of this field.
EXAMPLES:
sage: K.<a> = GF(9)
sage: K._cache.order()
9
Returns the order of this field.
EXAMPLES:
sage: K.<a> = GF(9)
sage: K._cache.order_c()
9
Return a random element of self.
EXAMPLES:
sage: k = GF(23**3, 'a')
sage: e = k._cache.random_element(); e
2*a^2 + 14*a + 21
sage: type(e)
<type 'sage.rings.finite_rings.element_givaro.FiniteField_givaroElement'>
sage: P.<x> = PowerSeriesRing(GF(3^3, 'a'))
sage: P.random_element(5)
2*a + 2 + (a^2 + a + 2)*x + (2*a + 1)*x^2 + (2*a^2 + a)*x^3 + 2*a^2*x^4 + O(x^5)
Bases: sage.rings.finite_rings.element_base.FinitePolyExtElement
An element of a (Givaro) finite field.
Return the string representation of self as an int (as returned by log_to_int()).
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: (b+1).int_repr()
'6'
Return the integer representation of self. When self is in the prime subfield, the integer returned is equal to self and not to log_repr.
Elements of this field are represented as ints in as follows:
for with
,
is
represented as:
.
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: k(4).integer_representation()
4
sage: b.integer_representation()
5
sage: type(b.integer_representation())
<type 'int'>
Return True if self == k(1).
EXAMPLES:
sage: k.<a> = GF(3^4); k
Finite Field in a of size 3^4
sage: a.is_one()
False
sage: k(1).is_one()
True
Return True if self is a square in self.parent()
ALGORITHM:
Elements are stored as powers of generators, so we simply check to see if it is an even power of a generator.
EXAMPLES:
sage: k.<a> = GF(9); k
Finite Field in a of size 3^2
sage: a.is_square()
False
sage: v = set([x^2 for x in k])
sage: [x.is_square() for x in v]
[True, True, True, True, True]
sage: [x.is_square() for x in k if not x in v]
[False, False, False, False]
TESTS:
sage: K = GF(27, 'a')
sage: set([a*a for a in K]) == set([a for a in K if a.is_square()])
True
sage: K = GF(25, 'a')
sage: set([a*a for a in K]) == set([a for a in K if a.is_square()])
True
sage: K = GF(16, 'a')
sage: set([a*a for a in K]) == set([a for a in K if a.is_square()])
True
Return True if self is nonzero, so it is a unit as an element of the finite field.
EXAMPLES:
sage: k.<a> = GF(3^4); k
Finite Field in a of size 3^4
sage: a.is_unit()
True
sage: k(0).is_unit()
False
Return the log to the base of self, i.e., an integer
such that
self.
Warning
TODO – This is currently implemented by solving the discrete log problem – which shouldn’t be needed because of how finite field elements are represented.
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: a = b^7
sage: a.log(b)
7
Return the log representation of self as a string. See the documentation of the _element_log_repr function of the parent field.
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: (b+2).log_repr()
'15'
Returns the int representation of self, as a Sage integer. Use int(self) to directly get a Python int.
Elements of this field are represented as ints in as follows:
for with
,
is
represented as:
.
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: k(4).log_to_int()
4
sage: b.log_to_int()
5
sage: type(b.log_to_int())
<type 'sage.rings.integer.Integer'>
Return the multiplicative order of this field element.
EXAMPLES:
sage: S.<b> = GF(5^2); S
Finite Field in b of size 5^2
sage: b.multiplicative_order()
24
sage: (b^6).multiplicative_order()
4
Return representation of this finite field element as a polynomial in the generator.
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: (b+2).poly_repr()
'b + 2'
Return self viewed as a polynomial over self.parent().prime_subfield().
EXAMPLES:
sage: k.<b> = GF(5^2); k
Finite Field in b of size 5^2
sage: f = (b^2+1).polynomial(); f
b + 4
sage: type(f)
<type 'sage.rings.polynomial.polynomial_zmod_flint.Polynomial_zmod_flint'>
sage: parent(f)
Univariate Polynomial Ring in b over Finite Field of size 5
Return a square root of this finite field element in its parent, if there is one. Otherwise, raise a ValueError.
INPUT:
extend – bool (default: True); if True, return a square root in an extension ring, if necessary. Otherwise, raise a ValueError if the root is not in the base ring.
Warning
this option is not implemented!
all – bool (default: False); if True, return all square roots of self, instead of just one.
Warning
The extend option is not implemented (yet).
ALGORITHM:
self is stored as for some generator
.
Return
for even
.
EXAMPLES:
sage: k.<a> = GF(7^2)
sage: k(2).sqrt()
3
sage: k(3).sqrt()
2*a + 6
sage: k(3).sqrt()**2
3
sage: k(4).sqrt()
2
sage: k.<a> = GF(7^3)
sage: k(3).sqrt()
Traceback (most recent call last):
...
ValueError: must be a perfect square.
TESTS:
sage: K = GF(49, 'a')
sage: all([a.sqrt()*a.sqrt() == a for a in K if a.is_square()])
True
sage: K = GF(27, 'a')
sage: all([a.sqrt()*a.sqrt() == a for a in K if a.is_square()])
True
sage: K = GF(8, 'a')
sage: all([a.sqrt()*a.sqrt() == a for a in K if a.is_square()])
True
sage: K.<a>=FiniteField(9)
sage: a.sqrt(extend = False, all = True)
[]
Bases: object
Iterator over FiniteField_givaro elements. We iterate multiplicatively, as powers of a fixed internal generator.
EXAMPLES:
sage: for x in GF(2^2,'a'): print x
0
a
a + 1
1
x.next() -> the next value, or raise StopIteration
EXAMPLES:
sage: k = GF(3**7, 'a')
sage: loads(dumps(k)) == k # indirect doctest
True
TESTS:
sage: k = GF(3**4, 'a')
sage: e = k.random_element()
sage: TestSuite(e).run() # indirect doctest