"Fossies" - the Fresh Open Source Software Archive

Member "NZMATH-1.2.0/nzmath/finitefield.py" (19 Nov 2012, 33292 Bytes) of package /linux/misc/old/NZMATH-1.2.0.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) Python source code syntax highlighting (style: standard) with prefixed line numbers. Alternatively you can here view or download the uninterpreted source code file. For more information about "finitefield.py" see the Fossies "Dox" file reference documentation.

    1 """
    2 finite fields.
    3 """
    4 
    5 from __future__ import division
    6 import logging
    7 import nzmath.gcd as gcd
    8 import nzmath.bigrandom as bigrandom
    9 import nzmath.bigrange as bigrange
   10 import nzmath.arith1 as arith1
   11 import nzmath.prime as prime
   12 import nzmath.ring as ring
   13 import nzmath.rational as rational
   14 import nzmath.factor.misc as factor_misc
   15 import nzmath.intresidue as intresidue
   16 import nzmath.poly.uniutil as uniutil
   17 import nzmath.matrix as matrix
   18 import nzmath.vector as vector
   19 import nzmath.compatibility
   20 
   21 _log = logging.getLogger('nzmath.finitefield')
   22 
   23 
   24 class FiniteFieldElement(ring.FieldElement):
   25     """
   26     The base class for all finite field element.
   27     """
   28     def __init__(self):
   29         # This class is abstract and can not be instantiated.
   30         if type(self) is FiniteFieldElement:
   31             raise NotImplementedError("the class FiniteFieldElement is abstract")
   32         ring.FieldElement.__init__(self)
   33 
   34 
   35 class FiniteField(ring.Field):
   36     """
   37     The base class for all finite fields.
   38     """
   39     def __init__(self, characteristic):
   40         # This class is abstract and can not be instantiated.
   41         if type(self) is FiniteField:
   42             raise NotImplementedError("the class FiniteField is abstract")
   43         ring.Field.__init__(self)
   44         self.char = characteristic
   45         self._orderfactor = None  # postpone the initialization
   46 
   47     def card(self):
   48         "Cardinality of the field"
   49         raise NotImplementedError("card have to be overridden")
   50 
   51     def getCharacteristic(self):
   52         """
   53         Return the characteristic of the field.
   54         """
   55         return self.char
   56 
   57     def order(self, elem):
   58         """
   59         Find and return the order of the element in the multiplicative
   60         group of the field.
   61         """
   62         if not elem:
   63             raise ValueError("zero is not in the multiplicative group.")
   64         if self._orderfactor is None:
   65             self._orderfactor = factor_misc.FactoredInteger(card(self) - 1)
   66         o = 1
   67         for p, e in self._orderfactor:
   68             b = elem ** (int(self._orderfactor) // (p**e))
   69             while b != self.one:
   70                 o = o * p
   71                 b = b ** p
   72         return o
   73 
   74     def random_element(self, *args):
   75         """
   76         Return a randomly chosen element og the field.
   77         """
   78         return self.createElement(bigrandom.randrange(*args))
   79 
   80     def primitive_element(self):
   81         """
   82         Return a primitive element of the field, i.e., a generator of
   83         the multiplicative group.
   84         """
   85         raise NotImplementedError("primitive_element should be overridden")
   86 
   87     def Legendre(self, element):
   88         """ Return generalize Legendre Symbol for FiniteField.
   89         """
   90         if not element:
   91             return 0
   92 
   93         if element == self.one or self.char == 2:
   94             return 1 # element of cyclic group order 2**n-1 also always 1
   95 
   96         # generic method:successive squaring
   97         # generalized Legendre symbol definition:
   98         #    (self/_ring) := self ** ((card(_ring)-1)/2)
   99         x = element ** ((card(self)-1) >> 1)
  100         if x == self.one:
  101             return 1
  102         elif x == -self.one:
  103             return -1
  104         raise ValueError("element must be not in field")
  105 
  106     def TonelliShanks(self, element):
  107         """ Return square root of element if exist.
  108         assume that characteristic have to be more than three.
  109         """
  110         if  self.char == 2:
  111             return self.sqrt(element) # should be error
  112 
  113         if self.Legendre(element) == -1:
  114             raise ValueError("There is no solution")
  115 
  116         # symbol and code reference from Cohen, CCANT 1.5.1
  117         (e, q) = arith1.vp(card(self)-1, 2)
  118 
  119         a = element
  120         n = self.createElement(self.char+1)
  121         while self.Legendre(n) != -1:
  122             n = self.random_element(2, card(self)) # field maybe large
  123         y = z = n ** q
  124         r = e
  125         x = a ** ((q-1) >> 1)
  126         b = a * (x ** 2)
  127         x = a * x
  128         while True:
  129             if b == self.one:
  130                 return x
  131             m = 1
  132             while m < r:
  133                 if b ** (2 ** m) == self.one:
  134                     break
  135                 m = m+1
  136             if m == r:
  137                 break
  138             t = y ** (2 ** (r-m-1))
  139             y = t ** 2
  140             r = m
  141             x = x * t
  142             b = b * y
  143         raise ValueError("There is no solution")
  144 
  145     def sqrt(self, element):
  146         """ Return square root if exist.
  147         """
  148         if not element or element == self.one:
  149             return element # trivial case
  150 
  151         # element of characteristic 2 always exist square root
  152         if  self.char == 2:
  153             return element ** ((card(self)) >> 1)
  154 
  155         # otherwise,
  156         return self.TonelliShanks(element)
  157 
  158 
  159 class FinitePrimeFieldElement(intresidue.IntegerResidueClass, FiniteFieldElement):
  160     """
  161     The class for finite prime field element.
  162     """
  163     def __init__(self, representative, modulus, modulus_is_prime=True):
  164         if not modulus_is_prime and not prime.primeq(abs(modulus)):
  165             raise ValueError("modulus must be a prime.")
  166 
  167         FiniteFieldElement.__init__(self)
  168         intresidue.IntegerResidueClass.__init__(self, representative, modulus)
  169 
  170         # ring
  171         self.ring = None
  172 
  173     def __repr__(self):
  174         return "FinitePrimeFieldElement(%d, %d)" % (self.n, self.m)
  175 
  176     def __str__(self):
  177         return "%d in F_%d" % (self.n, self.m)
  178 
  179     def getRing(self):
  180         """
  181         Return the finite prime field to which the element belongs.
  182         """
  183         if self.ring is None:
  184             self.ring = FinitePrimeField.getInstance(self.m)
  185         return self.ring
  186 
  187     def order(self):
  188         """
  189         Find and return the order of the element in the multiplicative
  190         group of F_p.
  191         """
  192         if self.n == 0:
  193             raise ValueError("zero is not in the group.")
  194         if not hasattr(self, "orderfactor"):
  195             self.orderfactor = factor_misc.FactoredInteger(self.m - 1)
  196         o = 1
  197         for p, e in self.orderfactor:
  198             b = self ** (int(self.orderfactor) // (p**e))
  199             while b.n != 1:
  200                 o = o*p
  201                 b = b**p
  202         return o
  203 
  204 
  205 class FinitePrimeField(FiniteField):
  206     """
  207     FinitePrimeField is also known as F_p or GF(p).
  208     """
  209 
  210     # class variable
  211     _instances = {}
  212 
  213     def __init__(self, characteristic):
  214         FiniteField.__init__(self, characteristic)
  215         self.registerModuleAction(rational.theIntegerRing, self._int_times)
  216         # mathematically Q_p = Q \ {r/s; gcd(r, s) == 1, gcd(s, p) > 1}
  217         # is more appropriate.
  218         self.registerModuleAction(rational.theRationalField, self._rat_times)
  219 
  220     @staticmethod
  221     def _int_times(integer, fpelem):
  222         """
  223         Return k * FinitePrimeFieldElement(n, p)
  224         """
  225         return FinitePrimeFieldElement(integer * fpelem.n, fpelem.m)
  226 
  227     @staticmethod
  228     def _rat_times(rat, fpelem):
  229         """
  230         Return Rational(a, b) * FinitePrimeFieldElement(n, p)
  231         """
  232         return FinitePrimeFieldElement(rat * fpelem.n, fpelem.m)
  233 
  234     def __eq__(self, other):
  235         if self is other:
  236             return True
  237         if isinstance(other, FinitePrimeField):
  238             return self.char == other.char
  239         return False
  240 
  241     def __hash__(self):
  242         return self.char
  243 
  244     def __ne__(self, other):
  245         return not (self == other)
  246 
  247     def __str__(self):
  248         return "F_%d" % self.char
  249 
  250     def __repr__(self):
  251         return "%s(%d)" % (self.__class__.__name__, self.char)
  252 
  253     def __hash__(self):
  254         return self.char & 0xFFFFFFFF
  255 
  256     def issubring(self, other):
  257         """
  258         Report whether another ring contains the field as a subring.
  259         """
  260         if self == other:
  261             return True
  262         if isinstance(other, FiniteField) and other.getCharacteristic() == self.char:
  263             return True
  264         try:
  265             return other.issuperring(self)
  266         except:
  267             return False
  268 
  269     def issuperring(self, other):
  270         """
  271         Report whether the field is a superring of another ring.
  272         Since the field is a prime field, it can be a superring of
  273         itself only.
  274         """
  275         return self == other
  276 
  277     def __contains__(self, elem):
  278         if isinstance(elem, FinitePrimeFieldElement) and elem.getModulus() == self.char:
  279             return True
  280         return False
  281 
  282     def createElement(self, seed):
  283         """
  284         Create an element of the field.
  285 
  286         'seed' should be an integer.
  287         """
  288         return FinitePrimeFieldElement(seed, self.char, modulus_is_prime=True)
  289 
  290     def primitive_element(self):
  291         """
  292         Return a primitive element of the field, i.e., a generator of
  293         the multiplicative group.
  294         """
  295         if self.char == 2:
  296             return self.one
  297         fullorder = card(self) - 1
  298         if self._orderfactor is None:
  299             self._orderfactor = factor_misc.FactoredInteger(fullorder)
  300         for i in bigrange.range(2, self.char):
  301             g = self.createElement(i)
  302             for p in self._orderfactor.prime_divisors():
  303                 if g ** (fullorder // p) == self.one:
  304                     break
  305             else:
  306                 return g
  307 
  308     def Legendre(self, element):
  309         """ Return generalize Legendre Symbol for FinitePrimeField.
  310         """
  311         if not element:
  312             return 0
  313         if element.n == 1 or element.m == 2:
  314             return 1 # trivial
  315 
  316         return arith1.legendre(element.n, element.m)
  317 
  318     def SquareRoot(self, element):
  319         """ Return square root if exist.
  320         """
  321         if not element or element.n == 1:
  322             return element # trivial case
  323         if element.m == 2:
  324             return element.getRing().one
  325         return arith1.modsqrt(element.n, element.m)
  326 
  327     def card(self):
  328         "Cardinality of the field"
  329         return self.char
  330 
  331     # properties
  332     def _getOne(self):
  333         "getter for one"
  334         if self._one is None:
  335             self._one = FinitePrimeFieldElement(1, self.char)
  336         return self._one
  337 
  338     one = property(_getOne, None, None, "multiplicative unit.")
  339 
  340     def _getZero(self):
  341         "getter for zero"
  342         if self._zero is None:
  343             self._zero = FinitePrimeFieldElement(0, self.char)
  344         return self._zero
  345 
  346     zero = property(_getZero, None, None, "additive unit.")
  347 
  348     @classmethod
  349     def getInstance(cls, characteristic):
  350         """
  351         Return an instance of the class with specified characteristic.
  352         """
  353         if characteristic not in cls._instances:
  354             cls._instances[characteristic] = cls(characteristic)
  355         return cls._instances[characteristic]
  356 
  357 
  358 FinitePrimeFieldPolynomial = uniutil.FinitePrimeFieldPolynomial
  359 uniutil.special_ring_table[FinitePrimeField] = FinitePrimeFieldPolynomial
  360 
  361 
  362 class ExtendedFieldElement(FiniteFieldElement):
  363     """
  364     ExtendedFieldElement is a class for an element of F_q.
  365     """
  366     def __init__(self, representative, field):
  367         """
  368         FiniteExtendedFieldElement(representative, field) creates
  369         an element of the finite extended field F_q.
  370 
  371         The argument representative has to be a polynomial with
  372         base field coefficients, i.e., if F_q is over F_{p^k}
  373         the representative has to F_{p^k} polynomial.
  374 
  375         Another argument field mut be an instance of
  376         ExtendedField.
  377         """
  378         if isinstance(field, ExtendedField):
  379             self.field = field
  380         else:
  381             raise TypeError("wrong type argument for field.")
  382         if (isinstance(representative, uniutil.FiniteFieldPolynomial) and
  383             representative.getCoefficientRing() == self.field.basefield):
  384             self.rep = self.field.modulus.mod(representative)
  385         else:
  386             _log.error(representative.__class__.__name__)
  387             raise TypeError("wrong type argument for representative.")
  388 
  389     def getRing(self):
  390         """
  391         Return the field to which the element belongs.
  392         """
  393         return self.field
  394 
  395     def _op(self, other, op):
  396         """
  397         Do `self (op) other'.
  398         op must be a name of the special method for binary operation.
  399         """
  400         if isinstance(other, ExtendedFieldElement):
  401             if other.field is self.field:
  402                 result = self.field.modulus.mod(getattr(self.rep, op)(other.rep))
  403                 return self.__class__(result, self.field)
  404             else:
  405                 fq1, fq2 = self.field, other.field
  406                 emb1, emb2 = double_embeddings(fq1, fq2)
  407                 return getattr(emb1(self), op)(emb2(other))
  408         if self.field.hasaction(ring.getRing(other)):
  409             # cases for action ring elements
  410             embedded = self.field.getaction(ring.getRing(other))(other, self.field.one)
  411             result = self.field.modulus.mod(getattr(self.rep, op)(embedded.rep))
  412             return self.__class__(result, self.field)
  413         else:
  414             return NotImplemented
  415 
  416     def __add__(self, other):
  417         """
  418         self + other
  419 
  420         other can be an element of either F_q, F_p or Z.
  421         """
  422         if other is 0 or other is self.field.zero:
  423             return self
  424         return self._op(other, "__add__")
  425 
  426     __radd__ = __add__
  427 
  428     def __sub__(self, other):
  429         """
  430         self - other
  431 
  432         other can be an element of either F_q, F_p or Z.
  433         """
  434         if other is 0 or other is self.field.zero:
  435             return self
  436         return self._op(other, "__sub__")
  437 
  438     def __rsub__(self, other):
  439         """
  440         other - self
  441 
  442         other can be an element of either F_q, F_p or Z.
  443         """
  444         return self._op(other, "__rsub__")
  445 
  446     def __mul__(self, other):
  447         """
  448         self * other
  449 
  450         other can be an element of either F_q, F_p or Z.
  451         """
  452         if other is 1 or other is self.field.one:
  453             return self
  454         return self._op(other, "__mul__")
  455 
  456     __rmul__ = __mul__
  457 
  458     def __truediv__(self, other):
  459         return self * other.inverse()
  460 
  461     __div__ = __truediv__
  462 
  463     def inverse(self):
  464         """
  465         Return the inverse of the element.
  466         """
  467         if not self:
  468             raise ZeroDivisionError("There is no inverse of zero.")
  469         return self.__class__(self.rep.extgcd(self.field.modulus)[0], self.field)
  470 
  471     def __pow__(self, index):
  472         """
  473         self ** index
  474 
  475         pow() with three arguments is not supported.
  476         """
  477         if not self:
  478             return self # trivial
  479         if index < 0:
  480             index %= self.field.getCharacteristic()
  481         power = pow(self.rep, index, self.field.modulus)
  482         return self.__class__(power, self.field)
  483 
  484     def __neg__(self):
  485         return self.field.zero - self
  486 
  487     def __pos__(self):
  488         return self
  489 
  490     def __eq__(self, other):
  491         try:
  492             if self.field == other.field:
  493                 if self.rep == other.rep:
  494                     return True
  495         except AttributeError:
  496             pass
  497         return False
  498 
  499     def __hash__(self):
  500         return hash(self.rep)
  501 
  502     def __ne__(self, other):
  503         return not (self == other)
  504 
  505     def __nonzero__(self):
  506         return bool(self.rep)
  507 
  508     def __repr__(self):
  509         return "%s(%s, %s)" % (self.__class__.__name__, repr(self.rep), repr(self.field))
  510 
  511     def trace(self):
  512         """
  513         Return the absolute trace.
  514         """
  515         p = self.field.char
  516         q = p
  517         alpha = self
  518         tr = alpha.rep[0]
  519         while q < card(self.field):
  520             q *= p
  521             alpha **= p
  522             tr += alpha.rep[0]
  523         if tr not in FinitePrimeField.getInstance(p):
  524             tr = FinitePrimeField.getInstance(p).createElement(tr)
  525         return tr
  526 
  527     def norm(self):
  528         """
  529         Return the absolute norm.
  530         """
  531         p = self.field.char
  532         nrm = (self ** ((card(self.field) - 1) // (p - 1))).rep[0]
  533         if nrm not in FinitePrimeField.getInstance(p):
  534             nrm = FinitePrimeField.getInstance(p).createElement(nrm)
  535         return nrm
  536 
  537 
  538 class ExtendedField(FiniteField):
  539     """
  540     ExtendedField is a class for finite field, whose cardinality
  541     q = p**n with a prime p and n>1. It is usually called F_q or GF(q).
  542     """
  543     def __init__(self, basefield, modulus):
  544         """
  545         ExtendedField(basefield, modulus)
  546 
  547         Create a field extension basefield[X]/(modulus(X)).
  548 
  549         The modulus has to be an irreducible polynomial with
  550         coefficients in the basefield.
  551         """
  552         FiniteField.__init__(self, basefield.char)
  553         self.basefield = basefield
  554         self.modulus = modulus
  555         if isinstance(self.basefield, FinitePrimeField):
  556             self.degree = self.modulus.degree()
  557         else:
  558             self.degree = self.basefield.degree * self.modulus.degree()
  559             bf = self.basefield
  560             while hasattr(bf, "basefield"):
  561                 # all basefields can be used as a scalar field.
  562                 self.registerModuleAction(bf.basefield, self._scalar_mul)
  563                 bf = bf.basefield
  564         self.registerModuleAction(self.basefield, self._scalar_mul)
  565         # integer is always a good scalar
  566         self.registerModuleAction(rational.theIntegerRing, self._scalar_mul)
  567 
  568     @classmethod
  569     def prime_extension(cls, characteristic, n_or_modulus):
  570         """
  571         ExtendedField.prime_extension(p, n_or_modulus) creates a
  572         finite field extended over prime field.
  573 
  574         characteristic must be prime. n_or_modulus can be:
  575           1) an integer greater than 1, or
  576           2) a polynomial in a polynomial ring of F_p with degree
  577              greater than 1.
  578         """
  579         if isinstance(n_or_modulus, (int, long)):
  580             n = n_or_modulus
  581             if n <= 1:
  582                 raise ValueError("degree of extension must be > 1.")
  583             # choose a method among three variants:
  584             #modulus = cls._random_irriducible(characteristic, n)
  585             #modulus = cls._small_irriducible(characteristic, n)
  586             modulus = cls._primitive_polynomial(characteristic, n)
  587         elif isinstance(n_or_modulus, FinitePrimeFieldPolynomial):
  588             modulus = n_or_modulus
  589             if not isinstance(modulus.getCoefficientRing(), FinitePrimeField):
  590                 raise TypeError("modulus must be F_p polynomial.")
  591             if modulus.degree() <= 1 or not modulus.isirreducible():
  592                 raise ValueError("modulus must be of degree greater than 1.")
  593         else:
  594             raise TypeError("degree or modulus must be supplied.")
  595         return cls(FinitePrimeField.getInstance(characteristic), modulus)
  596 
  597     @staticmethod
  598     def _random_irriducible(char, degree):
  599         """
  600         Return randomly chosen irreducible polynomial of self.degree.
  601         """
  602         cardinality = char ** degree
  603         basefield = FinitePrimeField.getInstance(char)
  604         seed = bigrandom.randrange(1, char) + cardinality
  605         cand = uniutil.FiniteFieldPolynomial(enumerate(arith1.expand(seed, cardinality)), coeffring=basefield)
  606         while cand.degree() < degree or not cand.isirreducible():
  607             seed = bigrandom.randrange(1, cardinality) + cardinality
  608             cand = uniutil.FiniteFieldPolynomial(enumerate(arith1.expand(seed, cardinality)), coeffring=basefield)
  609         _log.debug(cand.order.format(cand))
  610         return cand
  611 
  612     @staticmethod
  613     def _small_irriducible(char, degree):
  614         """
  615         Return an irreducible polynomial of self.degree with a small
  616         number of non-zero coefficients.
  617         """
  618         cardinality = char ** degree
  619         basefield = FinitePrimeField.getInstance(char)
  620         top = uniutil.polynomial({degree: 1}, coeffring=basefield)
  621         for seed in range(degree - 1):
  622             for const in range(1, char):
  623                 coeffs = [const] + arith1.expand(seed, 2)
  624                 cand = uniutil.polynomial(enumerate(coeffs), coeffring=basefield) + top
  625                 if cand.isirreducible():
  626                     _log.debug(cand.order.format(cand))
  627                     return cand
  628         for subdeg in range(degree):
  629             subseedbound = char ** subdeg
  630             for subseed in range(subseedbound + 1, char * subseedbound):
  631                 if not subseed % char:
  632                     continue
  633                 seed = subseed + cardinality
  634                 cand = uniutil.polynomial(enumerate(arith1.expand(seed, cardinality)), coeffring=basefield)
  635                 if cand.isirreducible():
  636                     return cand
  637 
  638     @staticmethod
  639     def _primitive_polynomial(char, degree):
  640         """
  641         Return a primitive polynomial of self.degree.
  642 
  643         REF: Lidl & Niederreiter, Introduction to finite fields and
  644              their applications.
  645         """
  646         cardinality = char ** degree
  647         basefield = FinitePrimeField.getInstance(char)
  648         const = basefield.primitive_element()
  649         if degree & 1:
  650             const = -const
  651         cand = uniutil.polynomial({0:const, degree:basefield.one}, basefield)
  652         maxorder = factor_misc.FactoredInteger((cardinality - 1) // (char - 1))
  653         var = uniutil.polynomial({1:basefield.one}, basefield)
  654         while not (cand.isirreducible() and
  655                    all(pow(var, int(maxorder) // p, cand).degree() > 0 for p in maxorder.prime_divisors())):
  656             # randomly modify the polynomial
  657             deg = bigrandom.randrange(1, degree)
  658             coeff = basefield.random_element(1, char)
  659             cand += uniutil.polynomial({deg:coeff}, basefield)
  660         _log.debug(cand.order.format(cand))
  661         return cand
  662 
  663     @staticmethod
  664     def _scalar_mul(integer, fqelem):
  665         """
  666         Return integer * (Fq element).
  667         """
  668         return fqelem.__class__(fqelem.rep * integer, fqelem.field)
  669 
  670     def card(self):
  671         """
  672         Return the cardinality of the field
  673         """
  674         return self.char ** self.degree
  675 
  676     def createElement(self, seed):
  677         """
  678         Create an element of the field.
  679         """
  680         if isinstance(seed, (int, long)):
  681             expansion = arith1.expand(seed, card(self.basefield))
  682             return ExtendedFieldElement(
  683                 uniutil.FiniteFieldPolynomial(enumerate(expansion), self.basefield),
  684                 self)
  685         elif seed in self.basefield:
  686             return ExtendedFieldElement(
  687                 uniutil.FiniteFieldPolynomial([(0, seed)], self.basefield),
  688                 self)
  689         elif seed in self:
  690             # seed is in self, return only embedding
  691             return self.zero + seed
  692         elif (isinstance(seed, uniutil.FiniteFieldPolynomial) and
  693               seed.getCoefficientRing() is self.basefield):
  694             return ExtendedFieldElement(seed, self)
  695         else:
  696             try:
  697                 # lastly check sequence
  698                 return ExtendedFieldElement(
  699                     uniutil.FiniteFieldPolynomial(enumerate(seed), self.basefield),
  700                     self)
  701             except TypeError:
  702                 raise TypeError("seed %s is not an appropriate object." % str(seed))
  703 
  704     def __repr__(self):
  705         return "%s(%d, %d)" % (self.__class__.__name__, self.char, self.degree)
  706 
  707     def __str__(self):
  708         return "F_%d @(%s)" % (card(self), str(self.modulus))
  709 
  710     def __hash__(self):
  711         return (self.char ** self.degree) & 0xFFFFFFFF
  712 
  713     def issuperring(self, other):
  714         """
  715         Report whether the field is a superring of another ring.
  716         """
  717         if self is other:
  718             return True
  719         elif isinstance(other, ExtendedField):
  720             if self.char == other.char and not (self.degree % other.degree):
  721                 return True
  722             return False
  723         elif isinstance(other, FinitePrimeField):
  724             if self.char == other.getCharacteristic():
  725                 return True
  726             return False
  727         try:
  728             return other.issubring(self)
  729         except:
  730             return False
  731 
  732     def issubring(self, other):
  733         """
  734         Report whether the field is a subring of another ring.
  735         """
  736         if self is other:
  737             return True
  738         elif isinstance(other, FinitePrimeField):
  739             return False
  740         elif isinstance(other, ExtendedField):
  741             if self.char == other.char and not (other.degree % self.degree):
  742                 return True
  743             return False
  744         try:
  745             return other.issuperring(self)
  746         except:
  747             return False
  748 
  749     def __contains__(self, elem):
  750         """
  751         Report whether elem is in field.
  752         """
  753         if (isinstance(elem, ExtendedFieldElement) and
  754             elem.getRing().modulus == self.modulus):
  755             return True
  756         elif (isinstance(elem, FinitePrimeFieldElement) and
  757               elem.getRing().getCharacteristic() == self.getCharacteristic()):
  758             return True
  759         return False
  760 
  761     def __eq__(self, other):
  762         """
  763         Equality test.
  764         """
  765         if isinstance(other, ExtendedField):
  766             return self.char == other.char and self.degree == other.degree
  767         return False
  768 
  769     def primitive_element(self):
  770         """
  771         Return a primitive element of the field, i.e., a generator of
  772         the multiplicative group.
  773         """
  774         fullorder = card(self) - 1
  775         if self._orderfactor is None:
  776             self._orderfactor = factor_misc.FactoredInteger(fullorder)
  777         for i in bigrange.range(self.char, card(self)):
  778             g = self.createElement(i)
  779             for p in self._orderfactor.prime_divisors():
  780                 if g ** (fullorder // p) == self.one:
  781                     break
  782             else:
  783                 return g
  784 
  785     # properties
  786     def _getOne(self):
  787         "getter for one"
  788         if self._one is None:
  789             self._one = ExtendedFieldElement(
  790                 uniutil.FiniteFieldPolynomial([(0, 1)], self.basefield),
  791                 self)
  792         return self._one
  793 
  794     one = property(_getOne, None, None, "multiplicative unit.")
  795 
  796     def _getZero(self):
  797         "getter for zero"
  798         if self._zero is None:
  799             self._zero = ExtendedFieldElement(
  800                 uniutil.FiniteFieldPolynomial([], self.basefield),
  801                 self)
  802         return self._zero
  803 
  804     zero = property(_getZero, None, None, "additive unit.")
  805 
  806 
  807 def fqiso(f_q, gfq):
  808     """
  809     Return isomorphism function of extended finite fields from f_q to gfq.
  810     """
  811     if f_q is gfq:
  812         return lambda x: x
  813     if card(f_q) != card(gfq):
  814         raise TypeError("both fields must have the same cardinality.")
  815 
  816     # find a root of f_q's defining polynomial in gfq.
  817     p = f_q.getCharacteristic()
  818     q = card(f_q)
  819     for i in bigrange.range(p, q):
  820         root = gfq.createElement(i)
  821         if not f_q.modulus(root):
  822             break
  823 
  824     # finally, define a function
  825     def f_q_to_gfq_iso(f_q_elem):
  826         """
  827         Return the image of the isomorphism of the given element.
  828         """
  829         if not f_q_elem:
  830             return gfq.zero
  831         if f_q_elem.rep.degree() == 0:
  832             # F_p elements
  833             return gfq.createElement(f_q_elem.rep)
  834         return f_q_elem.rep(root)
  835 
  836     return f_q_to_gfq_iso
  837 
  838 
  839 def embedding(f_q1, f_q2):
  840     """
  841     Return embedding homomorphism function from f_q1 to f_q2,
  842     where q1 = p ** k1, q2 = p ** k2 and k1 divides k2.
  843     """
  844     if card(f_q1) == card(f_q2):
  845         return fqiso(f_q1, f_q2)
  846     # search multiplicative generators of both fields and relate them.
  847     # 0. initialize basic variables
  848     q1, q2 = card(f_q1), card(f_q2)
  849 
  850     # 1. find a multiplicative generator of f_q2
  851     f_q2_gen = f_q2.primitive_element()
  852     f_q2_subgen = f_q2_gen ** ((q2 - 1) // (q1 - 1))
  853 
  854     # 2. find a root of defining polynomial of f_q1 in f_q2
  855     image_of_x_1 = _findroot(f_q1, f_q2, f_q2_subgen)
  856 
  857     # 3. finally, define a function
  858     def f_q1_to_f_q2_homo(f_q1_elem):
  859         """
  860         Return the image of the isomorphism of the given element.
  861         """
  862         if not f_q1_elem:
  863             return f_q2.zero
  864         if f_q1_elem.rep.degree() == 0:
  865             # F_p elements
  866             return f_q2.createElement(f_q1_elem.rep)
  867         return f_q1_elem.rep(image_of_x_1)
  868 
  869     return f_q1_to_f_q2_homo
  870 
  871 def double_embeddings(f_q1, f_q2):
  872     """
  873     Return embedding homomorphism functions from f_q1 and f_q2
  874     to the composite field.
  875     """
  876     identity = lambda x: x
  877     if f_q1 is f_q2:
  878         return (identity, identity)
  879     p = f_q2.getCharacteristic()
  880     k1, k2 = f_q1.degree, f_q2.degree
  881     if not k2 % k1:
  882         return (embedding(f_q1, f_q2), identity)
  883     if not k1 % k2:
  884         return (identity, embedding(f_q2, f_q1))
  885     composite = FiniteExtendedField(p, gcd.lcm(k1, k2))
  886     return (embedding(f_q1, composite), embedding(f_q2, composite))
  887 
  888 
  889 def _findroot(f_q1, f_q2, f_q2_subgen):
  890     """
  891     Find root of the defining polynomial of f_q1 in f_q2
  892     """
  893     if card(f_q1) > 50: # 50 is small, maybe
  894         _log.debug("by affine multiple (%d)" % card(f_q1))
  895         return affine_multiple_method(f_q1.modulus, f_q2)
  896     root = f_q2_subgen
  897     for i in range(1, card(f_q1)):
  898         if not f_q1.modulus(root):
  899             image_of_x_1 = root
  900             break
  901         root *= f_q2_subgen
  902     return image_of_x_1
  903 
  904 def affine_multiple_method(lhs, field):
  905     """
  906     Find and return a root of the equation lhs = 0 by brute force
  907     search in the given field.  If there is no root in the field,
  908     ValueError is raised.
  909 
  910     The first argument lhs is a univariate polynomial with
  911     coefficients in a finite field.  The second argument field is
  912     an extension field of the field of coefficients of lhs.
  913 
  914     Affine multiple A(X) is $\sum_{i=0}^{n} a_i X^{q^i} - a$ for some
  915     a_i's and a in the coefficient field of lhs, which is a multiple
  916     of the lhs.
  917     """
  918     polynomial_ring = lhs.getRing()
  919     coeff_field = lhs.getCoefficientRing()
  920     q = card(coeff_field)
  921     n = lhs.degree()
  922 
  923     # residues = [x, x^q, x^{q^2}, ..., x^{q^{n-1}}]
  924     residues = [lhs.mod(polynomial_ring.one.term_mul((1, 1)))] # x
  925     for i in range(1, n):
  926         residues.append(pow(residues[-1], q, lhs)) # x^{q^i}
  927 
  928     # find a linear relation among residues and a constant
  929     coeff_matrix = matrix.createMatrix(n, n, [coeff_field.zero] * (n**2), coeff_field)
  930     for j, residue in enumerate(residues):
  931         for i in range(residue.degree() + 1):
  932             coeff_matrix[i + 1, j + 1] = residue[i]
  933     constant_components = [coeff_field.one] + [coeff_field.zero] * (n - 1)
  934     constant_vector = vector.Vector(constant_components)
  935     try:
  936         relation_vector, kernel = coeff_matrix.solve(constant_vector)
  937         for j in range(n, 0, -1):
  938             if relation_vector[j]:
  939                 constant = relation_vector[j].inverse()
  940                 relation = [constant * c for c in relation_vector]
  941                 break
  942     except matrix.NoInverseImage:
  943         kernel_matrix = coeff_matrix.kernel()
  944         relation_vector = kernel_matrix[1]
  945         assert type(relation_vector) is vector.Vector
  946         for j in range(n, 0, -1):
  947             if relation_vector[j]:
  948                 normalizer = relation_vector[j].inverse()
  949                 relation = [normalizer * c for c in relation_vector]
  950                 constant = coeff_field.zero
  951                 break
  952 
  953     # define L(X) = A(X) + constant
  954     coeffs = {}
  955     for i, relation_i in enumerate(relation):
  956         coeffs[q**i] = relation_i
  957     linearized = uniutil.polynomial(coeffs, coeff_field)
  958 
  959     # Fq basis [1, X, ..., X^{s-1}]
  960     qbasis = [1]
  961     root = field.createElement(field.char)
  962     s = arith1.log(card(field), q)
  963     qbasis += [root**i for i in range(1, s)]
  964     # represent L as Matrix
  965     lmat = matrix.createMatrix(s, s, field.basefield)
  966     for j, base in enumerate(qbasis):
  967         imagei = linearized(base)
  968         if imagei.getRing() == field.basefield:
  969             lmat[1, j + 1] = imagei
  970         else:
  971             for i, coeff in imagei.rep.iterterms():
  972                 if coeff:
  973                     lmat[i + 1, j + 1] = coeff
  974     # solve L(X) = the constant
  975     constant_components = [constant] + [coeff_field.zero] * (s - 1)
  976     constant_vector = vector.Vector(constant_components)
  977     solution, kernel = lmat.solve(constant_vector)
  978     assert lmat * solution == constant_vector
  979     solutions = [solution]
  980     for v in kernel:
  981         for i in range(card(field.basefield)):
  982             solutions.append(solution + i * v)
  983 
  984     # roots of A(X) contains the solutions of lhs = 0
  985     for t in bigrange.multirange([(card(field.basefield),)] * len(kernel)):
  986         affine_root_vector = solution
  987         for i, ti in enumerate(t):
  988             affine_root_vector += ti * kernel[i]
  989         affine_root = field.zero
  990         for i, ai in enumerate(affine_root_vector):
  991             affine_root += ai * qbasis[i]
  992         if not lhs(affine_root):
  993             return affine_root
  994 
  995     raise ValueError("no root found")
  996 
  997 
  998 def FiniteExtendedField(characteristic, n_or_modulus):
  999     """
 1000     Return ExtendedField F_{p^n} or F_p[]/(modulus).
 1001 
 1002     This is a convenience wrapper for backward compatibility.
 1003     """
 1004     return ExtendedField.prime_extension(characteristic, n_or_modulus)
 1005 
 1006 FiniteExtendedFieldElement = ExtendedFieldElement