Ticket #14990: trac_14990-algebraic_closure_finite_field.patch

File trac_14990-algebraic_closure_finite_field.patch, 43.5 KB (added by pbruin, 6 years ago)

update (doctest coverage)

  • doc/en/reference/finite_rings/index.rst

    # HG changeset patch
    # User Peter Bruin <P.Bruin@warwick.ac.uk>
    # Date 1376472492 -7200
    # Node ID 7576e85707a44d974e5abe4827abedc97142a454
    # Parent  33a4ab61f20f40e2d3ba5f075ab97a4ac368b6b7
    Trac 14990: implement algebraic closures of finite fields
    
    diff --git a/doc/en/reference/finite_rings/index.rst b/doc/en/reference/finite_rings/index.rst
    a b  
    1515   sage/rings/finite_rings/finite_field_pari_ffelt
    1616   sage/rings/finite_rings/finite_field_prime_modn
    1717   sage/rings/finite_rings/homset
     18   sage/rings/algebraic_closure_finite_field
    1819
    1920
    2021.. include:: ../footer.txt
  • sage/matrix/matrix2.pyx

    diff --git a/sage/matrix/matrix2.pyx b/sage/matrix/matrix2.pyx
    a b  
    46964696            sage: A._eigenspace_format('galois') == 'galois'
    46974697            True
    46984698 
    4699         The algebraic closure of the rationals, the field of algebraic numbers,
    4700         (aka ``QQbar``) is implemented, while in general the algebraic closure
    4701         of a finite field is not implemented.  ::
     4699        The algebraic closure of the rationals (the field ``QQbar`` of
     4700        algebraic numbers) is implemented, as are algebraic closures
     4701        of finite fields::
    47024702       
    47034703            sage: A = matrix(QQ, 2, range(4))
    47044704            sage: A._eigenspace_format(None) == 'all'
    47054705            True
    47064706            sage: B = matrix(GF(13), 2, range(4))
    4707             sage: B._eigenspace_format(None) == 'galois'
    4708             True
     4707            sage: B._eigenspace_format(None)
     4708            'all'
    47094709       
    47104710        Subrings are promoted to fraction fields and then checked for the
    47114711        existence of algebraic closures.  ::
     
    47184718            msg = "format keyword must be None, 'all' or 'galois', not {0}"
    47194719            raise ValueError(msg.format(format))
    47204720
    4721         # For subrings of implemented algebraically closed fields we
    4722         #   default to all eigenspaces in the absence of a format keyword
    4723         # Add newly implemented algebraically closed fields to list below
    4724         #   and implement the determintion of the actual eigenvalues
    4725         #   in the eigenspace_left() routine
     4721        # In the absence of a format keyword, we default to 'all' for
     4722        # subrings of fields of which an algebraic closure is implemented.
    47264723        if format is None:
    4727             R = self.base_ring()
    4728             from sage.rings.qqbar import QQbar
    47294724            try:
    4730                 if R.fraction_field().algebraic_closure() in [QQbar]:
    4731                     return 'all'
     4725                F = self.base_ring().fraction_field()
     4726                if F.is_finite():
     4727                    _ = F.algebraic_closure('z')
    47324728                else:
    4733                     return 'galois'
     4729                    _ = F.algebraic_closure()
     4730                return 'all'
    47344731            except (NotImplementedError, AttributeError):
    47354732                return 'galois'
    47364733        else:
     
    1052110518            ValueError: matrix is not positive definite,
    1052210519            so cannot compute Cholesky decomposition           
    1052310520
    10524         Even in light of the above, you can sometimes get lucky
    10525         and arrive at a situation where a particular matrix has
    10526         a Cholesky decomposition when the general characteristics
    10527         of the matrix suggest this routine would fail. In this
    10528         example, the indefinite factorization produces a
    10529         diagonal matrix whose elements from the finite field
    10530         convert naturally to positive integers and are also
    10531         perfect squares.  ::
     10521        In certain cases, the algorithm can find an analogue of the
     10522        Cholesky decomposition over finite fields::
    1053210523
    1053310524            sage: F.<a> = FiniteField(5^3)
    1053410525            sage: A = matrix(F, [[         4,       2*a^2 + 3,         4*a + 1],
     
    1054010531            sage: L*L.transpose() == A
    1054110532            True
    1054210533
     10534            sage: F = FiniteField(7)
     10535            sage: A = matrix(F, [[4, 0], [0, 3]])
     10536            sage: A.cholesky()
     10537            [       2        0]
     10538            [       0 2*z2 + 6]
     10539
    1054310540        TESTS:
    1054410541
    1054510542        This verifies that :trac:`11274` is resolved.  ::
     
    1059010587                    sqrt = F(x.sqrt())
    1059110588                except (TypeError, ValueError):
    1059210589                    try:
    10593                         F = F.algebraic_closure()
     10590                        if F.is_finite():
     10591                            F = F.algebraic_closure('z')
     10592                        else:
     10593                            F = F.algebraic_closure()
    1059410594                    except (NotImplementedError, AttributeError):
    1059510595                        msg = "base field needs an algebraic closure with square roots, not {0}"
    1059610596                        raise TypeError(msg.format(F))
    10597                     # try again
    10598                     sqrt = F(x.sqrt())
     10597                    sqrt = F(x).sqrt()
    1059910598                splits.append(sqrt)
    1060010599            # move square root of the diagonal matrix
    1060110600            # into the lower triangular matrix
  • new file sage/rings/algebraic_closure_finite_field.py

    diff --git a/sage/rings/algebraic_closure_finite_field.py b/sage/rings/algebraic_closure_finite_field.py
    new file mode 100644
    - +  
     1r"""
     2Algebraic closures of finite fields
     3
     4Let `\Bold{F}` be a finite field, and let `\overline{\Bold{F}}` be an
     5algebraic closure of `\Bold{F}`; this is unique up to (non-canonical)
     6isomorphism.  For every `n\ge 1`, there is a unique subfield
     7`\Bold{F}_n` of `\overline{\Bold{F}}` such that
     8`\Bold{F}\subset\Bold{F}_n` and '[\Bold{F}_n:\Bold{F}]=1`.
     9
     10In Sage, algebraic closures of finite fields are implemented using
     11compatible systems of finite fields.  The resulting Sage object keeps
     12track of a finite lattice of the subfields `\Bold{F}_n` and the
     13embeddings between them.  This lattice is extended as necessary.
     14
     15The Sage class corresponding to `\overline{\Bold{F}}` can be
     16constructed from the finite field `\Bold{F}` by using the
     17:meth:`~sage.rings.finite_rings.finite_field_base.FiniteField.algebraic_closure`
     18method.  This invokes the ``AlgebraicClosureFiniteField`` factory
     19object to get unique representation.
     20
     21The Sage class for elements of `\overline{\Bold{F}}` is
     22:class:`AlgebraicClosureFiniteFieldElement`.  Such an element is
     23represented as an element of one of the `\Bold{F}_n`.  This means that
     24each element `x\in\Bold{F}` has infinitely many different
     25representations, one for each `n` such that `x` is in `\Bold{F}_n`.
     26
     27.. NOTE::
     28
     29    Only prime finite fields are currently accepted as base fields for
     30    algebraic closures.  To obtain an algebraic closure of a non-prime
     31    finite field `\Bold{F}`, take an algebraic closure of the prime
     32    field of `\Bold{F}` and embed `\Bold{F}` into this.
     33
     34    Algebraic closures of finite fields are currently implemented
     35    using (pseudo-)Conway polynomials; see
     36    :class:`AlgebraicClosureFiniteField_pseudo_conway` and the module
     37    :mod:`~sage.rings.finite_rings.conway_polynomials`.  Other
     38    implementations may be added by creating appropriate subclasses of
     39    :class:`AlgebraicClosureFiniteField_generic`.
     40
     41TEST::
     42
     43    sage: F = GF(5).algebraic_closure('z')
     44    sage: TestSuite(F).run()
     45
     46AUTHORS:
     47
     48- Peter Bruin (August 2013): initial version
     49
     50"""
     51
     52from sage.rings.finite_rings.element_base import is_FiniteFieldElement
     53from sage.rings.finite_rings.finite_field_base import is_FiniteField
     54from sage.rings.ring import Field
     55from sage.structure.element import FieldElement
     56from sage.structure.factory import UniqueFactory
     57
     58class AlgebraicClosureFiniteFieldElement(FieldElement):
     59    """
     60    Element of an algebraic closure of a finite field.
     61
     62    EXAMPLE::
     63
     64        sage: F = GF(3).algebraic_closure('z')
     65        sage: F.gen(2)
     66        z2
     67        sage: type(F.gen(2))
     68        <class 'sage.rings.algebraic_closure_finite_field.AlgebraicClosureFiniteField_pseudo_conway_with_category.element_class'>
     69
     70    """
     71    def __init__(self, parent, value):
     72        """
     73        TEST::
     74
     75            sage: F = GF(3).algebraic_closure('z')
     76            sage: TestSuite(F.gen(2)).run()
     77
     78        """
     79        if is_FiniteFieldElement(value):
     80            n = value.parent().degree()
     81        else:
     82            from sage.rings.integer import Integer
     83            n = Integer(1)
     84        self._value = parent._subfield(n).coerce(value)
     85        self._level = n
     86        FieldElement.__init__(self, parent)
     87
     88    def __reduce__(self):
     89        """
     90        Used for pickling.
     91
     92        TESTS::
     93
     94            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     95            sage: F = AlgebraicClosureFiniteField(GF(2), 'z')
     96            sage: loads(dumps(F.gen(5))) == F.gen(5)
     97            True
     98
     99        """
     100        return unpickle_AlgebraicClosureFiniteFieldElement, (self.parent(), self._level, self._repr_())
     101
     102    def _repr_(self):
     103        """
     104        Return a string representation of ``self``.
     105
     106        EXAMPLE::
     107
     108            sage: F = GF(3).algebraic_closure('z')
     109            sage: F._repr_()
     110            'Algebraic closure of Finite Field of size 3'
     111
     112        """
     113        return self._value._repr_()
     114
     115    def __cmp__(self, right):
     116        """
     117        Compare ``self`` with ``right``.
     118
     119        EXAMPLE::
     120
     121            sage: F = GF(3).algebraic_closure('z')
     122            sage: F.gen(2) == F.gen(3)
     123            False
     124
     125        """
     126        x, y = self.parent()._coerce_2(self, right)
     127        return cmp(x, y)
     128
     129    def _add_(self, right):
     130        """
     131        Return ``self`` + ``right``.
     132
     133        EXAMPLE::
     134
     135            sage: F = GF(3).algebraic_closure('z')
     136            sage: F.gen(2) + F.gen(3)
     137            z6^5 + 2*z6^4 + 2*z6^3 + z6^2 + 2*z6 + 1
     138
     139        """
     140        F = self.parent()
     141        x, y = F._coerce_2(self, right)
     142        return self.__class__(F, x + y)
     143
     144    def _sub_(self, right):
     145        """
     146        Return ``self`` - ``right``.
     147
     148        EXAMPLE::
     149
     150            sage: F = GF(3).algebraic_closure('z')
     151            sage: F.gen(2) - F.gen(3)
     152            z6^4 + 2*z6^3 + z6^2 + 2*z6
     153
     154        """
     155        F = self.parent()
     156        x, y = F._coerce_2(self, right)
     157        return self.__class__(F, x - y)
     158
     159    def _mul_(self, right):
     160        """
     161        Return ``self`` * ``right``.
     162
     163        EXAMPLE::
     164
     165            sage: F = GF(3).algebraic_closure('z')
     166            sage: F.gen(2) * F.gen(3)
     167            z6^5 + 2*z6^4 + z6^2 + 2
     168
     169        """
     170        F = self.parent()
     171        x, y = F._coerce_2(self, right)
     172        return self.__class__(F, x * y)
     173
     174    def _div_(self, right):
     175        """
     176        Return ``self`` / ``right``.
     177
     178        EXAMPLE::
     179
     180            sage: F = GF(3).algebraic_closure('z')
     181            sage: F.gen(2) / F.gen(3)
     182            z6^5 + 2*z6^4 + z6^3 + 1
     183
     184        """
     185        F = self.parent()
     186        x, y = F._coerce_2(self, right)
     187        return self.__class__(F, x / y)
     188
     189    def _change_level(self, n):
     190        """
     191        Return a representation of ``self`` as an element of the
     192        subfield of degree ``n`` of the parent, if possible.
     193
     194        EXAMPLE::
     195
     196            sage: F = GF(3).algebraic_closure('z')
     197            sage: z = F.gen(4)
     198            sage: (z^10)._change_level(6)
     199            2*z6^5 + 2*z6^3 + z6^2 + 2*z6 + 2
     200
     201        """
     202        F = self.parent()
     203        l = self._level
     204        m = l.gcd(n)
     205        xl = self._value
     206        xm = F.inclusion(m, l).section()(xl)
     207        xn = F.inclusion(m, n)(xm)
     208        return self.__class__(F, xn)
     209
     210    def is_square(self):
     211        """
     212        Return ``True`` if ``self`` is a square.
     213
     214        This always returns ``True``.
     215
     216        EXAMPLE::
     217
     218            sage: F = GF(3).algebraic_closure('z')
     219            sage: F.gen(2).is_square()
     220            True
     221
     222        """
     223        return True
     224
     225    def sqrt(self):
     226        """
     227        Return a square root of ``self``.
     228
     229        EXAMPLE::
     230
     231            sage: F = GF(3).algebraic_closure('z')
     232            sage: F.gen(2).sqrt()
     233            z4^3 + z4 + 1
     234
     235        """
     236        F = self.parent()
     237        x = self._value
     238        try:
     239            return self.__class__(F, x.sqrt(extend=False))
     240        except ValueError:
     241            l = self._level
     242            x = F.inclusion(l, 2*l)(x)
     243            return self.__class__(F, x.sqrt(extend=False))
     244
     245
     246def unpickle_AlgebraicClosureFiniteFieldElement(parent, level, x):
     247    """
     248    Unpickle an element `x` of an algebraic closure of a finite field.
     249
     250    TEST::
     251
     252        sage: F = GF(7).algebraic_closure('z')
     253        sage: loads(dumps(F.gen(2))) == F.gen(2)  # indirect doctest
     254        True
     255
     256    """
     257    return parent.coerce(parent._subfield(level)(x))
     258
     259
     260class AlgebraicClosureFiniteField_generic(Field):
     261    """
     262    Algebraic closure of a finite field.
     263    """
     264    def __init__(self, base_ring, name, category=None):
     265        """
     266        TEST::
     267
     268            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_generic
     269            sage: F = AlgebraicClosureFiniteField_generic(GF(5), 'z')
     270            sage: F
     271            Algebraic closure of Finite Field of size 5
     272
     273        """
     274        Field.__init__(self, base_ring=base_ring, names=name,
     275                       normalize=False, category=category)
     276
     277    def __getstate__(self):
     278        """
     279        Used for pickling.
     280
     281        TEST::
     282
     283            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_generic
     284            sage: F = AlgebraicClosureFiniteField_generic(GF(5), 'z')
     285            sage: F.__getstate__() is None
     286            True
     287
     288        """
     289        pass
     290
     291    def __setstate__(self, state):
     292        """
     293        Used for pickling.
     294
     295        TEST::
     296
     297            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_generic
     298            sage: F = AlgebraicClosureFiniteField_generic(GF(5), 'z')
     299            sage: F.__setstate__(None)
     300
     301        """
     302        pass
     303
     304    def __reduce__(self):
     305        """
     306        Used for pickling.
     307
     308        TEST::
     309
     310            sage: F = GF(5).algebraic_closure('z')
     311            sage: loads(dumps(F)) == F
     312            True
     313
     314        """
     315        return (AlgebraicClosureFiniteField,
     316                (self.base_ring(), self.variable_name(), self.category()),
     317                self.__getstate__())
     318
     319    def __cmp__(self, other):
     320        """
     321        Compare ``self`` with ``other``.
     322
     323        TEST::
     324
     325            sage: F3 = GF(3).algebraic_closure('z')
     326            sage: F3 == F3
     327            True
     328            sage: F5 = GF(5).algebraic_closure('z')
     329            sage: F3 == F5
     330            False
     331
     332        """
     333        if self is other:
     334            return 0
     335        c = cmp(type(self), type(other))
     336        if c != 0:
     337            return c
     338        return cmp((self.base_ring(), self.variable_name(), self.category()),
     339                   (other.base_ring(), other.variable_name(), other.category()))
     340
     341    def characteristic(self):
     342        """
     343        Return the characteristic of ``self``.
     344
     345        EXAMPLE::
     346
     347            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     348            sage: p = next_prime(1000)
     349            sage: F = AlgebraicClosureFiniteField(GF(p), 'z')
     350            sage: F.characteristic() == p
     351            True
     352
     353        """
     354        return self.base_ring().characteristic()
     355
     356    Element = AlgebraicClosureFiniteFieldElement
     357
     358    def _element_constructor_(self, x):
     359        """
     360        Construct an element of ``self``.
     361
     362        TEST::
     363
     364            sage: F = GF(5).algebraic_closure('z')
     365            sage: type(F(3))
     366            <class 'sage.rings.algebraic_closure_finite_field.AlgebraicClosureFiniteField_pseudo_conway_with_category.element_class'>
     367
     368        """
     369        if isinstance(x, self.element_class) and x.parent() is self:
     370            return x
     371        else:
     372            return self.element_class(self, x)
     373
     374    def _coerce_map_from_(self, other):
     375        """
     376        Return ``True`` if elements of ``other`` can be coerced into
     377        ``self``.
     378
     379        EXAMPLE::
     380
     381            sage: F = GF(7).algebraic_closure('z')
     382            sage: F.has_coerce_map_from(Integers())
     383            True
     384
     385        """
     386        if other is self:
     387            return True
     388        elif is_FiniteField(other) and self._subfield(other.degree()) is other:
     389            return True
     390        elif self._subfield(1).has_coerce_map_from(other):
     391            return True
     392
     393    def _repr_(self):
     394        """
     395        Return a string representation of ``self``.
     396
     397        EXAMPLE::
     398
     399            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     400            sage: F = AlgebraicClosureFiniteField(GF(5), 'z')
     401            sage: F._repr_()
     402            'Algebraic closure of Finite Field of size 5'
     403
     404        """
     405        return 'Algebraic closure of %s' % self.base_ring()
     406
     407    def _coerce_2(self, x, y):
     408        """
     409        Coerce `x` and `y` to a common subfield of ``self``.
     410
     411        TEST::
     412
     413            sage: F = GF(3).algebraic_closure('z')
     414            sage: x, y = F._coerce_2(F.gen(2), F.gen(3))
     415            sage: x.parent()
     416            Finite Field in z6 of size 3^6
     417            sage: y.parent()
     418            Finite Field in z6 of size 3^6
     419
     420        """
     421        n = x._level.lcm(y._level)
     422        mx = self.inclusion(x._level, n)
     423        my = self.inclusion(y._level, n)
     424        return mx(x._value), my(y._value)
     425
     426    def _get_polynomial(self, n):
     427        """
     428        Return the polynomial defining the unique subfield of degree
     429        `n` of ``self``.
     430
     431        This must be implemented by subclasses.
     432
     433        EXAMPLE::
     434
     435            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_generic
     436            sage: F = AlgebraicClosureFiniteField_generic(GF(5), 'z')
     437            sage: F._get_polynomial(1)
     438            Traceback (most recent call last):
     439            ...
     440            NotImplementedError
     441
     442        """
     443        raise NotImplementedError
     444
     445    def _get_im_gen(self, m, n):
     446        """
     447        Return the image of ``self.gen(m)`` under the canonical
     448        inclusion into ``self.subfield(n)``.
     449
     450        This must be implemented by subclasses.
     451
     452        EXAMPLE::
     453
     454            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_generic
     455            sage: F = AlgebraicClosureFiniteField_generic(GF(5), 'z')
     456            sage: F._get_im_gen(2, 4)
     457            Traceback (most recent call last):
     458            ...
     459            NotImplementedError
     460
     461        """
     462        raise NotImplementedError
     463
     464    def _subfield(self, n):
     465        """
     466        Return the unique subfield of degree `n` of ``self``.
     467
     468        EXAMPLE::
     469
     470            sage: F = GF(3).algebraic_closure('z')
     471            sage: F._subfield(4)
     472            Finite Field in z4 of size 3^4
     473
     474        """
     475        if n == 1:
     476            return self.base_ring()
     477        else:
     478            from sage.rings.finite_rings.constructor import FiniteField
     479            return FiniteField(self.base_ring().cardinality() ** n,
     480                               name=self.variable_name() + str(n),
     481                               modulus=self._get_polynomial(n))
     482
     483    def subfield(self, n):
     484        """
     485        Return the unique subfield of degree `n` of ``self``
     486        together with its canonical embedding into ``self``.
     487
     488        EXAMPLE::
     489
     490            sage: F = GF(3).algebraic_closure('z')
     491            sage: F.subfield(1)
     492            (Finite Field of size 3,
     493             Ring morphism:
     494               From: Finite Field of size 3
     495               To:   Algebraic closure of Finite Field of size 3
     496               Defn: 1 |--> 1)
     497            sage: F.subfield(4)
     498            (Finite Field in z4 of size 3^4,
     499             Ring morphism:
     500               From: Finite Field in z4 of size 3^4
     501               To:   Algebraic closure of Finite Field of size 3
     502               Defn: z4 |--> z4)
     503
     504        """
     505        Fn = self._subfield(n)
     506        return Fn, Fn.hom((self.gen(n),))
     507
     508    def inclusion(self, m, n):
     509        """
     510        Return the canonical inclusion map from the subfield
     511        of degree `m` to the subfield of degree `n`.
     512
     513        EXAMPLE::
     514
     515            sage: F = GF(3).algebraic_closure('z')
     516            sage: F.inclusion(1, 2)
     517            Ring Coercion morphism:
     518              From: Finite Field of size 3
     519              To:   Finite Field in z2 of size 3^2
     520            sage: F.inclusion(2, 4)
     521            Ring morphism:
     522              From: Finite Field in z2 of size 3^2
     523              To:   Finite Field in z4 of size 3^4
     524              Defn: z2 |--> 2*z4^3 + 2*z4^2 + 1
     525
     526        """
     527        if m == 1:
     528            return self.base_ring().hom(self._subfield(n))
     529        elif m.divides(n):
     530            return self._subfield(m).hom((self._get_im_gen(m, n),))
     531        else:
     532            raise ValueError("subfield of degree %s not contained in subfield of degree %s" % (m, n))
     533
     534    def ngens(self):
     535        """
     536        Return the number of generators of ``self``, which is
     537        infinity.
     538
     539        EXAMPLE::
     540
     541            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     542            sage: AlgebraicClosureFiniteField(GF(5), 'z').ngens()
     543            +Infinity
     544
     545        """
     546        from sage.rings.infinity import Infinity
     547        return Infinity
     548
     549    def gen(self, n):
     550        """
     551        Return the `n`-th generator of ``self``.
     552
     553        EXAMPLE::
     554
     555            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     556            sage: F = AlgebraicClosureFiniteField(GF(5), 'z')
     557            sage: F.gen(2)
     558            z2
     559
     560        """
     561        F = self._subfield(n)
     562        return self(F.gen())
     563
     564    def _first_ngens(self, n):
     565        """
     566        Return the first `n` generators of ``self``.
     567
     568        EXAMPLE::
     569
     570            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     571            sage: F = AlgebraicClosureFiniteField(GF(5), 'z')
     572            sage: F._first_ngens(3)
     573            (1, z2, z3)
     574
     575        """
     576        return tuple([self.gen(i + 1) for i in xrange(n)])
     577
     578    def algebraic_closure(self):
     579        """
     580        Return an algebraic closure of ``self``.
     581
     582        This always returns ``self``.
     583
     584        EXAMPLE::
     585
     586            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     587            sage: F = AlgebraicClosureFiniteField(GF(5), 'z')
     588            sage: F.algebraic_closure() is F
     589            True
     590
     591        """
     592        return self
     593
     594
     595class AlgebraicClosureFiniteField_pseudo_conway(AlgebraicClosureFiniteField_generic):
     596    """
     597    Algebraic closure of a finite field, constructed using
     598    pseudo-Conway polynomials.
     599
     600    EXAMPLE::
     601
     602        sage: F = GF(5).algebraic_closure('z')
     603        sage: type(F)
     604        <class 'sage.rings.algebraic_closure_finite_field.AlgebraicClosureFiniteField_pseudo_conway_with_category'>
     605
     606    """
     607    def __init__(self, base_ring, name, category=None):
     608        """
     609        TEST::
     610
     611            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_pseudo_conway
     612            sage: F = AlgebraicClosureFiniteField_pseudo_conway(GF(3), 'z')
     613            sage: F
     614            Algebraic closure of Finite Field of size 3
     615
     616        """
     617        if not (is_FiniteField(base_ring) and base_ring.is_prime_field()):
     618            raise NotImplementedError
     619        from sage.rings.finite_rings.conway_polynomials import PseudoConwayLattice
     620        self._pseudo_conway_lattice = PseudoConwayLattice(base_ring.characteristic())
     621        AlgebraicClosureFiniteField_generic.__init__(self, base_ring, name, category)
     622
     623    def __getstate__(self):
     624        """
     625        Used for pickling.
     626
     627        TEST::
     628
     629            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_pseudo_conway
     630            sage: F = AlgebraicClosureFiniteField_pseudo_conway(GF(5), 'z')
     631            sage: F.__getstate__()
     632            <class 'sage.rings.finite_rings.conway_polynomials.PseudoConwayLattice'>
     633
     634        """
     635        return self._pseudo_conway_lattice
     636
     637    def __setstate__(self, state):
     638        """
     639        Used for pickling.
     640
     641        TEST::
     642
     643            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_pseudo_conway
     644            sage: F = AlgebraicClosureFiniteField_pseudo_conway(GF(5), 'z')
     645            sage: F.__setstate__(F.__getstate__())
     646
     647        """
     648        self._pseudo_conway_lattice = state
     649
     650    def __cmp__(self, other):
     651        """
     652        Compare ``self`` with ``other``.
     653
     654        TEST::
     655
     656            sage: F3 = GF(3).algebraic_closure('z')
     657            sage: F3 == F3
     658            True
     659            sage: F5 = GF(5).algebraic_closure('z')
     660            sage: F3 == F5
     661            False
     662
     663        """
     664        c = AlgebraicClosureFiniteField_generic.__cmp__(self, other)
     665        if c != 0:
     666            return c
     667        return cmp(self._pseudo_conway_lattice, other._pseudo_conway_lattice)
     668
     669    def _get_polynomial(self, n):
     670        """
     671        Return the defining polynomial of the unique subfield of
     672        degree `n` of ``self``.
     673
     674        EXAMPLE::
     675
     676            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_pseudo_conway
     677            sage: F = AlgebraicClosureFiniteField_pseudo_conway(GF(5), 'z')
     678            sage: F._get_polynomial(1)
     679            x + 3
     680
     681        """
     682        return self._pseudo_conway_lattice.polynomial(n)
     683
     684    def _get_im_gen(self, m, n):
     685        """
     686        Return the image of ``self.gen(m)`` under the canonical
     687        inclusion into ``self.subfield(n)``.
     688
     689        EXAMPLE::
     690
     691            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField_pseudo_conway
     692            sage: F = AlgebraicClosureFiniteField_pseudo_conway(GF(5), 'z')
     693            sage: F._get_im_gen(2, 4)
     694            z4^3 + z4^2 + z4 + 3
     695
     696        """
     697        p = self.characteristic()
     698        return self._subfield(n).gen() ** ((p**n - 1)//(p**m - 1))
     699
     700
     701class AlgebraicClosureFiniteFieldFactory(UniqueFactory):
     702    """
     703    Factory for constructing algebraic closures of finite fields.
     704
     705    EXAMPLE::
     706
     707        sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     708        sage: F = GF(2).algebraic_closure('z')
     709        sage: F1 = AlgebraicClosureFiniteField(GF(2), 'z')
     710        sage: F1 is F
     711        True
     712        sage: loads(dumps(F)) is F
     713        True
     714
     715    """
     716    def create_key(self, base_ring, name, category=None, implementation=None, **kwds):
     717        """
     718        TEST::
     719
     720            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     721            sage: AlgebraicClosureFiniteField.create_key(GF(3), 'z')
     722            (Finite Field of size 3, 'z', Category of fields, 'pseudo_conway')
     723
     724        """
     725        if category is None:
     726            from sage.categories.fields import Fields
     727            category = Fields()
     728        if implementation is None:
     729            implementation = 'pseudo_conway'
     730        return (base_ring, name, category, implementation)
     731
     732    def create_object(self, version, key, **kwds):
     733        """
     734        TEST::
     735
     736            sage: from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     737            sage: key = (GF(3), 'z', Fields(), 'pseudo_conway')
     738            sage: AlgebraicClosureFiniteField.create_object(0, key)
     739            Algebraic closure of Finite Field of size 3
     740
     741        """
     742        base_ring, name, category, implementation = key
     743        if implementation == 'pseudo_conway':
     744            return AlgebraicClosureFiniteField_pseudo_conway(base_ring, name, category, **kwds)
     745        else:
     746            raise ValueError('unknown implementation for algebraic closure of finite field: %s'
     747                             % implementation)
     748
     749
     750AlgebraicClosureFiniteField = AlgebraicClosureFiniteFieldFactory('sage.rings.algebraic_closure_finite_field.AlgebraicClosureFiniteField')
  • sage/rings/finite_rings/constructor.py

    diff --git a/sage/rings/finite_rings/constructor.py b/sage/rings/finite_rings/constructor.py
    a b  
    261261        sage: K.<a> = GF(5**5, name='a', modulus=x^3 + 3*x + 3)
    262262        Traceback (most recent call last):
    263263        ...
    264         ValueError: The degree of the modulus does not correspond to the
    265         cardinality of the field.
     264        ValueError: the degree of the modulus does not equal the
     265        degree of the field.
    266266
    267267    If you wish to live dangerously, you can tell the constructor not
    268268    to test irreducibility using ``check_irreducible=False``, but this
     
    492492                    if not modulus.is_irreducible():
    493493                        raise ValueError("finite field modulus must be irreducible but it is not.")
    494494                    if modulus.degree() != n:
    495                         raise ValueError("The degree of the modulus does not correspond to the cardinality of the field.")
     495                        raise ValueError("the degree of the modulus does not equal the degree of the field.")
    496496                if name is None:
    497497                    raise TypeError("you must specify the generator name.")
    498498                if impl is None:
  • sage/rings/finite_rings/conway_polynomials.py

    diff --git a/sage/rings/finite_rings/conway_polynomials.py b/sage/rings/finite_rings/conway_polynomials.py
    a b  
    157157        else:
    158158            self.nodes = {}
    159159
     160    def __cmp__(self, other):
     161        """
     162        TEST::
     163
     164            sage: from sage.rings.finite_rings.conway_polynomials import PseudoConwayLattice
     165            sage: PCL3 = PseudoConwayLattice(3)
     166            sage: PCL5 = PseudoConwayLattice(5)
     167            sage: PCL3 == PCL3
     168            True
     169            sage: PCL3 == PCL5
     170            False
     171
     172        """
     173        if self is other:
     174            return 0
     175        c = cmp(type(self), type(other))
     176        if c != 0:
     177            return c
     178        return cmp(self.nodes, other.nodes)
     179
    160180    def polynomial(self, n):
    161181        r"""
    162182        Return the pseudo-Conway polynomial of degree `n` in this
  • sage/rings/finite_rings/element_pari_ffelt.pyx

    diff --git a/sage/rings/finite_rings/element_pari_ffelt.pyx b/sage/rings/finite_rings/element_pari_ffelt.pyx
    a b  
    193193                pari_catch_sig_on()
    194194                self.construct((<FiniteFieldElement_pari_ffelt>x).val)
    195195            else:
    196                 # This is where we *would* do coercion from one finite field to another...
    197196                raise TypeError("no coercion defined")
    198197
    199198        elif isinstance(x, Integer):
  • sage/rings/finite_rings/finite_field_base.pyx

    diff --git a/sage/rings/finite_rings/finite_field_base.pyx b/sage/rings/finite_rings/finite_field_base.pyx
    a b  
    1010"""
    1111include "sage/ext/stdsage.pxi"
    1212
     13from sage.categories.finite_fields import FiniteFields
    1314from sage.structure.parent cimport Parent
    1415from sage.misc.cachefunc import cached_method
    1516from sage.misc.prandom import randrange
     
    7475        """
    7576        return self
    7677
    77 from sage.categories.finite_fields import FiniteFields
    78 _FiniteFields = FiniteFields()
     78
    7979cdef class FiniteField(Field):
    8080    """
    8181    Abstract base class for finite fields.
    8282    """
    83     def __init__(self, base, names, normalize):
     83    def __init__(self, base, names, normalize, category=None):
    8484        """
    8585        Initialize ``self``.
    8686
     
    9797            sage: loads(K.dumps()) == K
    9898            True
    9999        """
    100         Field.__init__(self, base, names, normalize, category=_FiniteFields)
     100        if category is None:
     101            category = FiniteFields()
     102        Field.__init__(self, base, names, normalize, category)
    101103
    102104    def __repr__(self):
    103105        """
     
    320322            ...
    321323            TypeError: images do not define a valid homomorphism           
    322324        """
    323         if (self.characteristic() != codomain.characteristic()):
    324             raise ValueError, "no map from %s to %s"%(self, codomain)
    325         if (len(im_gens) != 1):
     325        if not codomain.characteristic().divides(self.characteristic()):
     326            raise ValueError, "no map from %s to %s" % (self, codomain)
     327        if len(im_gens) != 1:
    326328            raise ValueError, "only one generator for finite fields."
    327329
    328         return (im_gens[0].charpoly())(self.gen(0)).is_zero()
     330        return self.modulus()(im_gens[0]).is_zero()
    329331   
    330     def _Hom_(self, codomain, cat=None):
     332    def _Hom_(self, codomain, category=None):
    331333        """
    332         Return homset of homomorphisms from ``self`` to the finite field
    333         codomain. This function is implicitly called by the Hom method or
     334        Return the set of homomorphisms from ``self`` to ``codomain``
     335        in ``category``.
     336
     337        This function is implicitly called by the ``Hom`` method or
    334338        function.
    335339
    336         The ``cat`` option is currently ignored.
    337 
    338340        EXAMPLES::
    339341
    340342            sage: K.<a> = GF(25); K
     
    343345            Automorphism group of Finite Field in a of size 5^2
    344346        """
    345347        from sage.rings.finite_rings.homset import FiniteFieldHomset
    346         return FiniteFieldHomset(self, codomain)
     348        from sage.rings.homset import RingHomset
     349        if category.is_subcategory(FiniteFields()):
     350            return FiniteFieldHomset(self, codomain, category)
     351        else:
     352            return RingHomset(self, codomain, category)
    347353
    348354    def gen(self):
    349355        r"""
     
    522528            sage: GF(997).order()
    523529            997
    524530        """
    525         raise NotImplementedError
     531        return self.characteristic()**self.degree()
    526532
    527533    def factored_order(self):
    528534        """
     
    804810            if R.characteristic() == self.characteristic():
    805811                if R.degree() == 1:
    806812                    return True
    807                 if self.degree() % R.degree() == 0:
    808                     if hasattr(self, '_prefix') and hasattr(R, '_prefix'):
    809                         return R.hom((self.gen() ** ((self.order() - 1)//(R.order() - 1)),))
     813                elif (R.degree().divides(self.degree())
     814                      and hasattr(self, '_prefix') and hasattr(R, '_prefix')):
     815                    return R.hom((self.gen() ** ((self.order() - 1)//(R.order() - 1)),))
    810816
    811817    def construction(self):
    812818        """
     
    979985                raise ValueError, "name must be None, a string or a dictionary indexed by divisors of the degree"
    980986            return [self.subfields(m, name=name[m])[0] for m in divisors]
    981987
    982     def algebraic_closure(self):
     988    def algebraic_closure(self, name):
    983989        """
    984         Return the algebraic closure of ``self`` (not implemented).
     990        Return an algebraic closure of ``self``.
     991
     992        INPUT:
     993
     994        - ``name`` -- string: prefix to use for variable names of
     995          subfields
     996
     997        EXAMPLE::
     998
     999            sage: F = GF(5).algebraic_closure('z')
     1000            sage: F
     1001            Algebraic closure of Finite Field of size 5
     1002            sage: F.gen(3)
     1003            z3
    9851004
    9861005        .. NOTE::
    9871006
    988            This is not yet implemented for finite fields.
     1007            This is currently only implemented for prime fields.
    9891008
    990         EXAMPLES::
    991 
    992             sage: GF(5).algebraic_closure()
    993             Traceback (most recent call last):
    994             ...
    995             NotImplementedError: Algebraic closures of finite fields not implemented.
    9961009        """
    997         raise NotImplementedError, "Algebraic closures of finite fields not implemented."
     1010        from sage.rings.algebraic_closure_finite_field import AlgebraicClosureFiniteField
     1011        return AlgebraicClosureFiniteField(self, name)
    9981012
    9991013    @cached_method
    10001014    def is_conway(self):
  • sage/rings/finite_rings/finite_field_ext_pari.py

    diff --git a/sage/rings/finite_rings/finite_field_ext_pari.py b/sage/rings/finite_rings/finite_field_ext_pari.py
    a b  
    520520                # canonically isomorphic finite fields
    521521                return element_ext_pari.FiniteField_ext_pariElement(self, x)
    522522            else:
    523                 # This is where we *would* do coercion from one finite field to another...
    524523                raise TypeError, "no coercion defined"
    525524
    526525        elif sage.interfaces.gap.is_GapElement(x):
  • sage/rings/finite_rings/finite_field_pari_ffelt.py

    diff --git a/sage/rings/finite_rings/finite_field_pari_ffelt.py b/sage/rings/finite_rings/finite_field_pari_ffelt.py
    a b  
    138138
    139139        self._modulus = modulus
    140140        self._degree = n
    141         self._card = p ** n
    142141        self._kwargs = {}
    143142
    144143        self._gen_pari = pari(modulus).ffgen()
     
    162161        try:
    163162            return self.__hash
    164163        except AttributeError:
    165             self.__hash = hash((self._card, self.variable_name(), self._modulus))
     164            self.__hash = hash((self.cardinality(), self.variable_name(), self._modulus))
    166165            return self.__hash
    167166
    168167    def __reduce__(self):
     
    197196        """
    198197        if not isinstance(other, FiniteField_pari_ffelt):
    199198            return cmp(type(self), type(other))
    200         return cmp((self._card, self.variable_name(), self._modulus),
    201                    (other._card, other.variable_name(), other._modulus))
     199        return cmp((self.cardinality(), self.variable_name(), self._modulus),
     200                   (other.cardinality(), other.variable_name(), other._modulus))
    202201
    203202    def __richcmp__(left, right, op):
    204203        """
     
    499498            return x
    500499        else:
    501500            return self.element_class(self, x)
    502 
    503     def order(self):
    504         """
    505         The number of elements of the finite field.
    506 
    507         EXAMPLE::
    508 
    509             sage: k = FiniteField(2^10, 'a', impl='pari_ffelt')
    510             sage: k
    511             Finite Field in a of size 2^10
    512             sage: k.order()
    513             1024
    514         """
    515         return self._card
  • sage/rings/finite_rings/finite_field_prime_modn.py

    diff --git a/sage/rings/finite_rings/finite_field_prime_modn.py b/sage/rings/finite_rings/finite_field_prime_modn.py
    a b  
    3535_FiniteFields = FiniteFields()
    3636
    3737import sage.rings.finite_rings.integer_mod_ring as integer_mod_ring
    38 import sage.rings.integer as integer
     38from sage.rings.integer import Integer
    3939import sage.rings.finite_rings.integer_mod as integer_mod
    40 import sage.rings.arith as arith
    4140
    4241from sage.rings.integer_ring import ZZ
    4342from sage.rings.finite_rings.integer_mod_ring import IntegerModRing_generic
     
    5655        sage: FiniteField(next_prime(1000))
    5756        Finite Field of size 1009
    5857    """
    59     def __init__(self, p, name=None, check=True):
     58    def __init__(self, p, check=True):
    6059        """
    6160        Return a new finite field of order `p` where `p` is prime.
    6261
     
    6463
    6564        - ``p`` -- an integer at least 2
    6665
    67         - ``name`` -- ignored
    68 
    6966        - ``check`` -- bool (default: ``True``); if ``False``, do not
    7067          check ``p`` for primality
    7168
     
    7471            sage: F = FiniteField(3); F
    7572            Finite Field of size 3
    7673        """
    77         p = integer.Integer(p)
    78         if check and not arith.is_prime(p):
     74        p = Integer(p)
     75        if check and not p.is_prime():
    7976            raise ArithmeticError, "p must be prime"
    8077        self.__char = p
    81         self._IntegerModRing_generic__factored_order = factorization.Factorization([(p,1)], integer.Integer(1))
     78        self._IntegerModRing_generic__factored_order = factorization.Factorization([(p,1)], Integer(1))
    8279        self._kwargs = {}
    8380        # FiniteField_generic does nothing more than IntegerModRing_generic, and
    8481        # it saves a non trivial overhead
     
    10198        r"""
    10299        Compare ``self`` with ``other``.
    103100
    104         Two finite prime fields are considered equal if their characteristic
    105         is equal.
     101        Two finite prime fields are considered equal if and only if
     102        their characteristics are equal.
    106103
    107104        EXAMPLES::
    108            
     105
    109106            sage: K = FiniteField(3)
    110107            sage: copy(K) == K
    111108            True
    112109        """
    113110        if not isinstance(other, FiniteField_prime_modn):
    114111            return cmp(type(self), type(other))
    115 #        elif other.__class__ != FiniteField_prime_modn:
    116 #            return -cmp(other, self)
    117112        return cmp(self.__char, other.__char)
    118113
    119114    def __richcmp__(left, right, op):
     
    214209
    215210    def modulus(self):
    216211        """
    217         Return the minimal polynomial of ``self``, which is always `x - 1`.
     212        Return the defining polynomial of ``self``.
     213
     214        This always returns `x - 1`.
    218215
    219216        EXAMPLES::
    220217
     
    284281
    285282    def gen(self, n=0):
    286283        """
    287         Return generator of this finite field as an extension of its
    288         prime field.
     284        Return a generator of ``self`` over its prime field.
     285
     286        This always returns 1.
    289287
    290288        .. NOTE::
    291289
     
    304302        """
    305303        if n != 0:
    306304            raise IndexError, "only one generator"
    307         return self(1)
     305        return self.one_element()
    308306
    309307    def __iter__(self):
    310308        """
     
    332330        while i:
    333331            yield i
    334332            i += one
    335        
     333
    336334    def degree(self):
    337335        """
    338         Returns the degree of the finite field, which is a positive
    339         integer.
     336        Return the degree of ``self`` over its prime field.
     337
     338        This always returns 1.
    340339
    341340        EXAMPLES::
    342341
    343342            sage: FiniteField(3).degree()
    344343            1
    345344        """
    346         return integer.Integer(1)
     345        return Integer(1)
  • sage/rings/finite_rings/integer_mod.pyx

    diff --git a/sage/rings/finite_rings/integer_mod.pyx b/sage/rings/finite_rings/integer_mod.pyx
    a b  
    349349        """
    350350        return sage.rings.finite_rings.integer_mod.mod, (self.lift(), self.modulus(), self.parent())
    351351       
     352    def _im_gens_(self, codomain, im_gens):
     353        """
     354        Return the image of ``self`` under the map that sends the
     355        generators of the parent to ``im_gens``.
     356
     357        EXAMPLE::
     358
     359            sage: a = Mod(7, 10)
     360            sage: R = ZZ.quotient(5)
     361            sage: a._im_gens_(R, (R(1),))
     362            2
     363        """
     364        return codomain._coerce_(self)
     365
    352366    def is_nilpotent(self):
    353367        r"""
    354368        Return ``True`` if ``self`` is nilpotent,
  • sage/rings/residue_field.pyx

    diff --git a/sage/rings/residue_field.pyx b/sage/rings/residue_field.pyx
    a b  
    13301330            #<class 'sage.rings.residue_field.ResidueFiniteField_prime_modn_with_category'>
    13311331        """
    13321332        ResidueField_generic.__init__(self, p)
    1333         FiniteField_prime_modn.__init__(self, intp, name)
     1333        FiniteField_prime_modn.__init__(self, intp)
    13341334        from sage.rings.finite_rings.integer_mod import IntegerMod_to_IntegerMod, Integer_to_IntegerMod, Int_to_IntegerMod
    13351335        K = OK = p.ring()
    13361336        if OK.is_field():
  • sage/structure/element.pyx

    diff --git a/sage/structure/element.pyx b/sage/structure/element.pyx
    a b  
    19261926    return IS_INSTANCE(x, CommutativeRingElement)
    19271927
    19281928cdef class CommutativeRingElement(RingElement):
    1929     def _im_gens_(self, codomain, im_gens):
    1930         if len(im_gens) == 1 and self._parent.gen(0) == 1:
    1931             return codomain(self)
    1932         raise NotImplementedError
    1933    
     1929    """
     1930    Base class for elements of commutative rings.
     1931    """
    19341932    def inverse_mod(self, I):
    19351933        r"""
    19361934        Return an inverse of self modulo the ideal `I`, if defined,