source: sage/rings/number_field/number_field_element.pyx @ 6387:a00907b28f6e

Revision 6387:a00907b28f6e, 65.5 KB checked in by Robert Bradshaw <robertwb@…>, 6 years ago (diff)

NumberFieldElement?_quadratic work
Line 
1"""
2Number Field Elements
3
4AUTHORS:
5    -- William Stein version before it got cython'd
6    -- Joel B. Mohler (2007-03-09): First reimplementation into cython
7    -- William Stein (2007-09-04): add doctests
8    -- Robert Bradshaw (2007-09-15): specialized classes for relative and absolute elements
9"""
10
11# TODO -- relative extensions need to be completely rewritten, so one
12# can get easy access to representation of elements in their relative
13# form.  Functions like matrix below can't be done until relative
14# extensions are re-written this way.  Also there needs to be class
15# hierarchy for number field elements, integers, etc.  This is a
16# nontrivial project, and it needs somebody to attack it.  I'm amazed
17# how long this has gone unattacked.
18
19# Relative elements need to be a derived class or something.  This is
20# terrible as it is now.
21
22#*****************************************************************************
23#       Copyright (C) 2004, 2007 William Stein <wstein@gmail.com>
24#
25#  Distributed under the terms of the GNU General Public License (GPL)
26#
27#    This code is distributed in the hope that it will be useful,
28#    but WITHOUT ANY WARRANTY; without even the implied warranty of
29#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
30#    General Public License for more details.
31#
32#  The full text of the GPL is available at:
33#
34#                  http://www.gnu.org/licenses/
35#*****************************************************************************
36
37import operator
38
39include "../../ext/stdsage.pxi"
40cdef extern from *:
41    # TODO: move to stdsage.pxi
42    object PY_NEW_SAME_TYPE(object o)
43    bint PY_TYPE_CHECK_EXACT(object o, object type)
44include '../../ext/interrupt.pxi'
45include '../../ext/python_int.pxi'
46
47import sage.rings.field_element
48import sage.rings.infinity
49import sage.rings.polynomial.polynomial_element
50import sage.rings.polynomial.polynomial_ring
51import sage.rings.rational_field
52import sage.rings.rational
53import sage.rings.integer_ring
54import sage.rings.integer
55import sage.rings.arith
56
57import sage.rings.number_field.number_field
58
59from sage.libs.ntl.ntl cimport ntl_ZZ, ntl_ZZX
60from sage.rings.integer_ring cimport IntegerRing_class
61
62from sage.libs.all import pari_gen
63from sage.libs.pari.gen import PariError
64
65QQ = sage.rings.rational_field.QQ
66ZZ = sage.rings.integer_ring.ZZ
67Integer_sage = sage.rings.integer.Integer
68
69from sage.rings.integer cimport Integer
70from sage.rings.rational cimport Rational
71
72
73def is_NumberFieldElement(x):
74    """
75    Return True if x is of type NumberFieldElement, i.e., an
76    element of a number field.
77
78    EXAMPLES:
79        sage: is_NumberFieldElement(2)
80        False
81        sage: k.<a> = NumberField(x^7 + 17*x + 1)
82        sage: is_NumberFieldElement(a+1)
83        True
84    """
85    return PY_TYPE_CHECK(x, NumberFieldElement)
86
87def __create__NumberFieldElement_version0(parent, poly):
88    """
89    Used in unpickling elements of number fields.
90
91    EXAMPLES:
92    Since this is just used in unpickling, we unpickle.
93   
94        sage: k.<a> = NumberField(x^3 - 2)
95        sage: loads(dumps(a+1)) == a + 1
96        True       
97    """
98    return NumberFieldElement(parent, poly)
99
100def __create__NumberFieldElement_version1(parent, cls, poly):
101    """
102    Used in unpickling elements of number fields.
103
104    EXAMPLES:
105    Since this is just used in unpickling, we unpickle.
106   
107        sage: k.<a> = NumberField(x^3 - 2)
108        sage: loads(dumps(a+1)) == a + 1
109        True       
110    """
111    return cls(parent, poly)
112
113cdef class NumberFieldElement(FieldElement):
114    """
115    An element of a number field.
116
117    EXAMPLES:
118        sage: k.<a> = NumberField(x^3 + x + 1)
119        sage: a^3
120        -a - 1       
121    """
122    cdef NumberFieldElement _new(self):
123        """
124        Quickly creates a new initialized NumberFieldElement with the
125        same parent as self.
126        """
127        cdef NumberFieldElement x
128        x = <NumberFieldElement>PY_NEW_SAME_TYPE(self)
129        x._parent = self._parent
130        return x
131
132    def __init__(self, parent, f):
133        """
134        INPUT:
135            parent -- a number field
136            f -- defines an element of a number field.
137
138        EXAMPLES:
139        The following examples illustrate creation of elements of
140        number fields, and some basic arithmetic.
141
142        First we define a polynomial over Q.
143            sage: R.<x> = PolynomialRing(QQ)
144            sage: f = x^2 + 1
145
146        Next we use f to define the number field.
147            sage: K.<a> = NumberField(f); K
148            Number Field in a with defining polynomial x^2 + 1
149            sage: a = K.gen()
150            sage: a^2
151            -1
152            sage: (a+1)^2
153            2*a
154            sage: a^2
155            -1
156            sage: z = K(5); 1/z
157            1/5
158
159        We create a cube root of 2.
160            sage: K.<b> = NumberField(x^3 - 2)
161            sage: b = K.gen()
162            sage: b^3
163            2
164            sage: (b^2 + b + 1)^3
165            12*b^2 + 15*b + 19
166
167        This example illustrates save and load:
168            sage: K.<a> = NumberField(x^17 - 2)
169            sage: s = a^15 - 19*a + 3
170            sage: loads(s.dumps()) == s
171            True
172        """
173        sage.rings.field_element.FieldElement.__init__(self, parent)
174
175        cdef ZZ_c coeff
176        if isinstance(f, (int, long, Integer_sage)):
177            # set it up and exit immediately
178            # fast pathway
179            (<Integer>ZZ(f))._to_ZZ(&coeff)
180            SetCoeff( self.__numerator, 0, coeff )
181            conv_ZZ_int( self.__denominator, 1 )
182            return
183           
184        elif isinstance(f, NumberFieldElement):
185            self.__numerator = (<NumberFieldElement>f).__numerator
186            self.__denominator = (<NumberFieldElement>f).__denominator
187            return
188
189        ppr = parent.polynomial_ring()
190        if isinstance(parent, sage.rings.number_field.number_field.NumberField_relative):
191            ppr = parent.base_field().polynomial_ring()
192
193        if isinstance(f, pari_gen):
194            f = f.lift()
195            f = ppr(f)
196        if not isinstance(f, sage.rings.polynomial.polynomial_element.Polynomial):
197            f = ppr(f)
198        if f.degree() >= parent.absolute_degree():
199            if f.variable_name() != 'x':
200                f = f.change_variable_name('x')
201            if isinstance(parent, sage.rings.number_field.number_field.NumberField_relative):
202                f %= parent.absolute_polynomial()
203            else:
204                f %= parent.polynomial()
205
206        # Set Denominator
207        den = f.denominator()
208        (<Integer>ZZ(den))._to_ZZ(&self.__denominator)
209
210        cdef long i
211        num = f * den
212        for i from 0 <= i <= num.degree():
213            (<Integer>ZZ(num[i]))._to_ZZ(&coeff)
214            SetCoeff( self.__numerator, i, coeff )
215
216    def __alloc__(self):
217        ZZX_construct(&self.__numerator)
218        ZZ_construct(&self.__denominator)
219
220    def __dealloc__(self):
221        ZZX_destruct(&self.__numerator)
222        ZZ_destruct(&self.__denominator)
223
224    def _lift_cyclotomic_element(self, new_parent):
225        """
226        Creates an element of the passed field from this field.  This
227        is specific to creating elements in a cyclotomic field from
228        elements in another cyclotomic field.  This function aims to
229        make this common coercion extremely fast!
230
231        EXAMPLES:
232            sage: C.<zeta5>=CyclotomicField(5)
233            sage: CyclotomicField(10)(zeta5+1)  # The function _lift_cyclotomic_element does the heavy lifting in the background
234            zeta10^2 + 1
235            sage: (zeta5+1)._lift_cyclotomic_element(CyclotomicField(10))  # There is rarely a purpose to call this function directly
236            zeta10^2 + 1
237
238        AUTHOR:
239            Joel B. Mohler
240        """
241        # Right now, I'm a little confused why quadratic extension fields have a zeta_order function
242        # I would rather they not have this function since I don't want to do this isinstance check here.
243        if not isinstance(self.parent(), sage.rings.number_field.number_field.NumberField_cyclotomic) or not isinstance(new_parent, sage.rings.number_field.number_field.NumberField_cyclotomic):
244            raise TypeError, "The field and the new parent field must both be cyclotomic fields."
245
246        try:
247            small_order = self.parent().zeta_order()
248            large_order = new_parent.zeta_order()
249        except AttributeError:
250            raise TypeError, "The field and the new parent field must both be cyclotomic fields."
251
252        try:
253            _rel = ZZ(large_order / small_order)
254        except TypeError:
255            raise TypeError, "The zeta_order of the new field must be a multiple of the zeta_order of the original."
256
257        cdef NumberFieldElement x = <NumberFieldElement>PY_NEW_SAME_TYPE(self)
258        x._parent = <ParentWithBase>new_parent
259        x.__denominator = self.__denominator
260        cdef ZZX_c result
261        cdef ZZ_c tmp
262        cdef int i
263        cdef int rel = _rel
264        cdef ntl_ZZX _num
265        cdef ntl_ZZ _den
266        _num, _den = new_parent.polynomial_ntl()
267        for i from 0 <= i <= deg(self.__numerator):
268            GetCoeff(tmp, self.__numerator, i)
269            SetCoeff(result, i*rel, tmp)
270        rem_ZZX(x.__numerator, result, _num.x[0])
271        return x
272
273    def __reduce__(self):
274        """
275        Used in pickling number field elements.
276
277        EXAMPLES:
278            sage: k.<a> = NumberField(x^3 - 17*x^2 + 1)
279            sage: t = a.__reduce__(); t
280            (<built-in function __create__NumberFieldElement_version1>, (Number Field in a with defining polynomial x^3 - 17*x^2 + 1, <type 'sage.rings.number_field.number_field_element.NumberFieldElement_absolute'>, x))
281            sage: t[0](*t[1]) == a
282            True       
283        """
284        return __create__NumberFieldElement_version1, \
285               (self.parent(), type(self), self.polynomial())
286
287    def __repr__(self):
288        """
289        String representation of this number field element,
290        which is just a polynomial in the generator.
291
292        EXAMPLES:
293            sage: k.<a> = NumberField(x^2 + 2)
294            sage: b = (2/3)*a + 3/5
295            sage: b.__repr__()
296            '2/3*a + 3/5'       
297        """
298        x = self.polynomial()
299        return str(x).replace(x.parent().variable_name(),self.parent().variable_name())
300
301    def _im_gens_(self, codomain, im_gens):
302        """
303        This is used in computing homomorphisms between number fields.
304
305        EXAMPLES:
306            sage: k.<a> = NumberField(x^2 - 2)
307            sage: m.<b> = NumberField(x^4 - 2)
308            sage: phi = k.hom([b^2])
309            sage: phi(a+1)
310            b^2 + 1
311            sage: (a+1)._im_gens_(m, [b^2])
312            b^2 + 1       
313        """
314        # NOTE -- if you ever want to change this so relative number fields are
315        # in terms of a root of a poly.
316        # The issue is that elements of a relative number field are represented in terms
317        # of a generator for the absolute field.  However the morphism gives the image
318        # of gen, which need not be a generator for the absolute field.  The morphism
319        # has to be *over* the relative element.
320        return codomain(self.polynomial()(im_gens[0]))
321
322    def _latex_(self):
323        """
324        Returns the latex representation for this element.
325
326        EXAMPLES:
327            sage: C,zeta12=CyclotomicField(12).objgen()
328            sage: latex(zeta12^4-zeta12)
329            \zeta_{12}^{2} - \zeta_{12} - 1
330        """
331        return self.polynomial()._latex_(name=self.parent().latex_variable_name())
332
333    def _pari_(self, var='x'):
334        raise NotImplementedError, "NumberFieldElement sub-classes must override _pari_"
335
336    def _pari_init_(self, var='x'):
337        """
338        Return GP/PARI string representation of self. This is used for
339        converting this number field element to GP/PARI.  The returned
340        string defines a pari Mod in the variable is var, which is by
341        default 'x' -- not the name of the generator of the number
342        field.
343
344        INPUT:
345            var -- (default: 'x') the variable of the pari Mod.
346
347        EXAMPLES:
348            sage: K.<a> = NumberField(x^5 - x - 1)
349            sage: ((1 + 1/3*a)^4)._pari_init_()
350            'Mod(1/81*x^4 + 4/27*x^3 + 2/3*x^2 + 4/3*x + 1, x^5 - x - 1)'
351            sage: ((1 + 1/3*a)^4)._pari_init_('a')
352            'Mod(1/81*a^4 + 4/27*a^3 + 2/3*a^2 + 4/3*a + 1, a^5 - a - 1)'
353
354        Note that _pari_init_ can fail because of reserved words in PARI,
355        and since it actually works by obtaining the PARI representation
356        of something.
357            sage: K.<theta> = NumberField(x^5 - x - 1)
358            sage: b = (1/2 - 2/3*theta)^3; b
359            -8/27*theta^3 + 2/3*theta^2 - 1/2*theta + 1/8
360            sage: b._pari_init_('theta')
361            Traceback (most recent call last):
362            ...
363            PariError: unexpected character (2)
364
365        Fortunately pari_init returns everything in terms of x by default.
366            sage: pari(b)
367            Mod(-8/27*x^3 + 2/3*x^2 - 1/2*x + 1/8, x^5 - x - 1)
368        """
369        return repr(self._pari_(var=var))
370##         if var == None:
371##             var = self.parent().variable_name()
372##         if isinstance(self.parent(), sage.rings.number_field.number_field.NumberField_relative):
373##             f = self.polynomial()._pari_()
374##             g = str(self.parent().pari_relative_polynomial())
375##             base = self.parent().base_ring()
376##             gsub = base.gen()._pari_()
377##             gsub = str(gsub).replace(base.variable_name(), "y")
378##             g = g.replace("y", gsub)
379##         else:
380##             f = str(self.polynomial()).replace("x",var)
381##             g = str(self.parent().polynomial()).replace("x",var)
382##         return 'Mod(%s, %s)'%(f,g)
383
384    def __getitem__(self, n):
385        """
386        Return the n-th coefficient of this number field element, written
387        as a polynomial in the generator.
388
389        Note that $n$ must be between 0 and $d-1$, where $d$ is the
390        degree of the number field.
391
392        EXAMPLES:
393            sage: m.<b> = NumberField(x^4 - 1789)
394            sage: c = (2/3-4/5*b)^3; c
395            -64/125*b^3 + 32/25*b^2 - 16/15*b + 8/27
396            sage: c[0]
397            8/27
398            sage: c[2]
399            32/25
400            sage: c[3]
401            -64/125
402
403        We illustrate bounds checking:
404            sage: c[-1]
405            Traceback (most recent call last):
406            ...
407            IndexError: index must be between 0 and degree minus 1.
408            sage: c[4]
409            Traceback (most recent call last):
410            ...
411            IndexError: index must be between 0 and degree minus 1.
412
413        The list method implicitly calls __getitem__:
414            sage: list(c)
415            [8/27, -16/15, 32/25, -64/125]
416            sage: m(list(c)) == c
417            True           
418        """
419        if n < 0 or n >= self.parent().degree():     # make this faster.
420            raise IndexError, "index must be between 0 and degree minus 1."
421        return self.polynomial()[n]
422
423    cdef int _cmp_c_impl(left, sage.structure.element.Element right) except -2:
424        cdef NumberFieldElement _right = right
425        return not (ZZX_equal(&left.__numerator, &_right.__numerator) and ZZ_equal(&left.__denominator, &_right.__denominator))
426
427    def __abs__(self):
428        r"""
429        Return the numerical absolute value of this number field
430        element with respect to the first archimedean embedding, to 53
431        bits of precision.
432
433        This is the \code{abs( )} Python function.  If you want a different
434        embedding or precision, use \code{self.abs(...)}.
435
436        EXAMPLES:
437            sage: k.<a> = NumberField(x^3 - 2)
438            sage: abs(a)
439            1.25992104989487
440            sage: abs(a)^3
441            2.00000000000000
442            sage: a.abs(prec=128)
443            1.2599210498948731647672106072782283506       
444        """
445        return self.abs(prec=53, i=0)
446
447    def abs(self, prec=53, i=0):
448        """
449        Return the absolute value of this element with respect to the
450        ith complex embedding of parent, to the given precision.
451
452        INPUT:
453            prec -- (default: 53) integer bits of precision
454            i -- (default: ) integer, which embedding to use
455
456        EXAMPLES:
457            sage: z = CyclotomicField(7).gen()
458            sage: abs(z)
459            1.00000000000000
460            sage: abs(z^2 + 17*z - 3)
461            16.0604426799931
462            sage: K.<a> = NumberField(x^3+17)
463            sage: abs(a)
464            2.57128159065824
465            sage: a.abs(prec=100)
466            2.5712815906582353554531872087
467            sage: a.abs(prec=100,i=1)
468            2.5712815906582353554531872087
469            sage: a.abs(100, 2)
470            2.5712815906582353554531872087
471
472        Here's one where the absolute value depends on the embedding.
473            sage: K.<b> = NumberField(x^2-2)
474            sage: a = 1 + b
475            sage: a.abs(i=0)
476            2.41421356237309
477            sage: a.abs(i=1)
478            0.414213562373095
479        """
480        P = self.parent().complex_embeddings(prec)[i]
481        return abs(P(self))
482
483    def coordinates_in_terms_of_powers(self):
484        r"""
485        Let $\alpha$ be self.  Return a Python function that takes any
486        element of the parent of self in $\QQ(\alpha)$ and writes it in
487        terms of the powers of $\alpha$: $1, \alpha, \alpha^2, ...$.
488
489        (NOT CACHED).
490
491        EXAMPLES:
492        This function allows us to write elements of a number field in
493        terms of a different generator without having to construct a
494        whole separate number field.
495
496            sage: y = polygen(QQ,'y'); K.<beta> = NumberField(y^3 - 2); K
497            Number Field in beta with defining polynomial y^3 - 2
498            sage: alpha = beta^2 + beta + 1
499            sage: c = alpha.coordinates_in_terms_of_powers(); c
500            Coordinate function that writes elements in terms of the powers of beta^2 + beta + 1
501            sage: c(beta)
502            [-2, -3, 1]
503            sage: c(alpha)
504            [0, 1, 0]
505            sage: c((1+beta)^5)
506            [3, 3, 3]
507            sage: c((1+beta)^10)
508            [54, 162, 189]
509       
510        This function works even if self only generates a subfield
511        of this number field.
512            sage: k.<a> = NumberField(x^6 - 5)
513            sage: alpha = a^3
514            sage: c = alpha.coordinates_in_terms_of_powers()
515            sage: c((2/3)*a^3 - 5/3)
516            [-5/3, 2/3]
517            sage: c
518            Coordinate function that writes elements in terms of the powers of a^3
519            sage: c(a)
520            Traceback (most recent call last):
521            ...
522            ArithmeticError: vector is not in free module
523        """
524        K = self.parent()
525        V, from_V, to_V = K.absolute_vector_space()
526        h = K(1)
527        B = [to_V(h)]
528        f = self.minpoly()
529        for i in range(f.degree()-1):
530            h *= self
531            B.append(to_V(h))
532        W = V.span_of_basis(B)
533        return CoordinateFunction(self, W, to_V)
534
535    def complex_embeddings(self, prec=53):
536        """
537        Return the images of this element in the floating point
538        complex numbers, to the given bits of precision.
539
540        INPUT:
541            prec -- integer (default: 53) bits of precision
542
543        EXAMPLES:
544            sage: k.<a> = NumberField(x^3 - 2)
545            sage: a.complex_embeddings()
546            [1.25992104989487, -0.629960524947437 + 1.09112363597172*I, -0.629960524947437 - 1.09112363597172*I]
547            sage: a.complex_embeddings(10)
548            [1.3, -0.63 + 1.1*I, -0.63 - 1.1*I]
549            sage: a.complex_embeddings(100)
550            [1.2599210498948731647672106073, -0.62996052494743658238360530364 + 1.0911236359717214035600726142*I, -0.62996052494743658238360530364 - 1.0911236359717214035600726142*I]
551        """
552        phi = self.parent().complex_embeddings(prec)
553        return [f(self) for f in phi]
554
555    def complex_embedding(self, prec=53, i=0):
556        """
557        Return the i-th embedding of self in the complex numbers, to
558        the given precision.
559
560        EXAMPLES:
561            sage: k.<a> = NumberField(x^3 - 2)
562            sage: a.complex_embedding()
563            1.25992104989487
564            sage: a.complex_embedding(10)
565            1.3
566            sage: a.complex_embedding(100)
567            1.2599210498948731647672106073
568            sage: a.complex_embedding(20, 1)
569            -0.62996 + 1.0911*I
570            sage: a.complex_embedding(20, 2)
571            -0.62996 - 1.0911*I
572        """
573        return self.parent().complex_embeddings(prec)[i](self)
574
575    def is_square(self, root=False):
576        """
577        Return True if self is a square in its parent number field and
578        otherwise return False.
579
580        INPUT:
581            root -- if True, also return a square root (or None if self
582                    is not a perfect square)
583
584        EXAMPLES:
585            sage: m.<b> = NumberField(x^4 - 1789)
586            sage: b.is_square()
587            False
588            sage: c = (2/3*b + 5)^2; c
589            4/9*b^2 + 20/3*b + 25
590            sage: c.is_square()
591            True
592            sage: c.is_square(True)
593            (True, 2/3*b + 5)
594
595        We also test the functional notation.
596            sage: is_square(c, True)
597            (True, 2/3*b + 5)
598            sage: is_square(c)
599            True
600            sage: is_square(c+1)
601            False
602        """
603        v = self.sqrt(all=True)
604        t = len(v) > 0
605        if root:
606            if t:
607                return t, v[0]
608            else:
609                return False, None
610        else:
611            return t
612       
613    def sqrt(self, all=False):
614        """
615        Returns the square root of this number in the given number field.
616       
617        EXAMPLES:
618            sage: K.<a> = NumberField(x^2 - 3)
619            sage: K(3).sqrt()
620            a
621            sage: K(3).sqrt(all=True)
622            [a, -a]
623            sage: K(a^10).sqrt()
624            9*a
625            sage: K(49).sqrt()
626            7
627            sage: K(1+a).sqrt()
628            Traceback (most recent call last):
629            ...
630            ValueError: a + 1 not a square in Number Field in a with defining polynomial x^2 - 3
631            sage: K(0).sqrt()
632            0
633            sage: K((7+a)^2).sqrt(all=True)
634            [a + 7, -a - 7]
635           
636            sage: K.<a> = CyclotomicField(7)
637            sage: a.sqrt() 
638            a^4
639           
640            sage: K.<a> = NumberField(x^5 - x + 1)
641            sage: (a^4 + a^2 - 3*a + 2).sqrt()
642            a^3 - a^2
643
644        ALGORITHM:
645            Use Pari to factor $x^2$ - \code{self} in K.
646           
647        """
648        # For now, use pari's factoring abilities
649        R = sage.rings.polynomial.polynomial_ring.PolynomialRing(self._parent, 't')
650        f = R([-self, 0, 1])
651        roots = f.roots()
652        if all:
653            return [r[0] for r in roots]
654        elif len(roots) > 0:
655            return roots[0][0]
656        else:
657            raise ValueError, "%s not a square in %s"%(self, self._parent)
658
659    cdef void _reduce_c_(self):
660        """
661        Pull out common factors from the numerator and denominator!
662        """
663        cdef ZZ_c gcd
664        cdef ZZ_c t1
665        cdef ZZX_c t2
666        content(t1, self.__numerator)
667        GCD_ZZ(gcd, t1, self.__denominator)
668        if sign(gcd) != sign(self.__denominator):
669            negate(t1, gcd)
670            gcd = t1
671        div_ZZX_ZZ(t2, self.__numerator, gcd)
672        div_ZZ_ZZ(t1, self.__denominator, gcd)
673        self.__numerator = t2
674        self.__denominator = t1
675
676    cdef ModuleElement _add_c_impl(self, ModuleElement right):
677        cdef NumberFieldElement x
678        cdef NumberFieldElement _right = right
679        x = self._new()
680        mul_ZZ(x.__denominator, self.__denominator, _right.__denominator)
681        cdef ZZX_c t1, t2
682        mul_ZZX_ZZ(t1, self.__numerator, _right.__denominator)
683        mul_ZZX_ZZ(t2, _right.__numerator, self.__denominator)
684        add_ZZX(x.__numerator, t1, t2)
685        x._reduce_c_()
686        return x
687
688    cdef ModuleElement _sub_c_impl(self, ModuleElement right):
689        cdef NumberFieldElement x
690        cdef NumberFieldElement _right = right
691        x = self._new()
692        mul_ZZ(x.__denominator, self.__denominator, _right.__denominator)
693        cdef ZZX_c t1, t2
694        mul_ZZX_ZZ(t1, self.__numerator, _right.__denominator)
695        mul_ZZX_ZZ(t2, _right.__numerator, self.__denominator)
696        sub_ZZX(x.__numerator, t1, t2)
697        x._reduce_c_()
698        return x
699
700    cdef RingElement _mul_c_impl(self, RingElement right):
701        """
702        Returns the product of self and other as elements of a number field.
703
704        EXAMPLES:
705            sage: C.<zeta12>=CyclotomicField(12)
706            sage: zeta12*zeta12^11
707            1
708            sage: G.<a> = NumberField(x^3 + 2/3*x + 1)
709            sage: a^3
710            -2/3*a - 1
711            sage: a^3+a
712            1/3*a - 1
713        """
714        cdef NumberFieldElement x
715        cdef NumberFieldElement _right = right
716        cdef ZZX_c temp
717        cdef ZZ_c temp1
718        cdef ZZ_c parent_den
719        cdef ZZX_c parent_num
720        self._parent_poly_c_( &parent_num, &parent_den )
721        x = self._new()
722        _sig_on
723        # MulMod doesn't handle non-monic polynomials.
724        # Therefore, we handle the non-monic case entirely separately.
725        if ZZX_is_monic( &parent_num ):
726            mul_ZZ(x.__denominator, self.__denominator, _right.__denominator)
727            MulMod_ZZX(x.__numerator, self.__numerator, _right.__numerator, parent_num)
728        else:
729            mul_ZZ(x.__denominator, self.__denominator, _right.__denominator)
730            mul_ZZX(x.__numerator, self.__numerator, _right.__numerator)
731            if ZZX_degree(&x.__numerator) >= ZZX_degree(&parent_num):
732                mul_ZZX_ZZ( x.__numerator, x.__numerator, parent_den )
733                mul_ZZX_ZZ( temp, parent_num, x.__denominator )
734                power_ZZ(temp1,LeadCoeff_ZZX(temp),ZZX_degree(&x.__numerator)-ZZX_degree(&parent_num)+1)
735                PseudoRem_ZZX(x.__numerator, x.__numerator, temp)
736                mul_ZZ(x.__denominator, x.__denominator, parent_den)
737                mul_ZZ(x.__denominator, x.__denominator, temp1)
738        _sig_off
739        x._reduce_c_()
740        return x
741
742        #NOTES: In LiDIA, they build a multiplication table for the
743        #number field, so it's not necessary to reduce modulo the
744        #defining polynomial every time:
745        #     src/number_fields/algebraic_num/order.cc: compute_table
746        # but asymptotically fast poly multiplication means it's
747        # actually faster to *not* build a table!?!
748
749    cdef RingElement _div_c_impl(self, RingElement right):
750        """
751        Returns the quotient of self and other as elements of a number field.
752
753        EXAMPLES:
754            sage: C.<I>=CyclotomicField(4)
755            sage: 1/I
756            -I
757            sage: I/0
758            Traceback (most recent call last):
759            ...
760            ZeroDivisionError: Number field element division by zero
761           
762            sage: G.<a> = NumberField(x^3 + 2/3*x + 1)
763            sage: a/a
764            1
765            sage: 1/a
766            -a^2 - 2/3
767            sage: a/0
768            Traceback (most recent call last):
769            ...
770            ZeroDivisionError: Number field element division by zero
771        """
772        cdef NumberFieldElement x
773        cdef NumberFieldElement _right = right
774        cdef ZZX_c inv_num
775        cdef ZZ_c inv_den
776        cdef ZZ_c parent_den
777        cdef ZZX_c parent_num
778        cdef ZZX_c temp
779        cdef ZZ_c temp1
780        if not _right:
781            raise ZeroDivisionError, "Number field element division by zero"
782        self._parent_poly_c_( &parent_num, &parent_den )
783        x = self._new()
784        _sig_on
785        _right._invert_c_(&inv_num, &inv_den)
786        if ZZX_is_monic( &parent_num ):
787            mul_ZZ(x.__denominator, self.__denominator, inv_den)
788            MulMod_ZZX(x.__numerator, self.__numerator, inv_num, parent_num)
789        else:
790            mul_ZZ(x.__denominator, self.__denominator, inv_den)
791            mul_ZZX(x.__numerator, self.__numerator, inv_num)
792            if ZZX_degree(&x.__numerator) >= ZZX_degree(&parent_num):
793                mul_ZZX_ZZ( x.__numerator, x.__numerator, parent_den )
794                mul_ZZX_ZZ( temp, parent_num, x.__denominator )
795                power_ZZ(temp1,LeadCoeff_ZZX(temp),ZZX_degree(&x.__numerator)-ZZX_degree(&parent_num)+1)
796                PseudoRem_ZZX(x.__numerator, x.__numerator, temp)
797                mul_ZZ(x.__denominator, x.__denominator, parent_den)
798                mul_ZZ(x.__denominator, x.__denominator, temp1)
799        x._reduce_c_()
800        _sig_off
801        return x
802
803    def __floordiv__(self, other):
804        """
805        Return the quotient of self and other.  Since these are field
806        elements the floor division is exactly the same as usual
807        division.
808
809        EXAMPLES:
810            sage: m.<b> = NumberField(x^4 + x^2 + 2/3)
811            sage: c = (1+b) // (1-b); c
812            3/4*b^3 + 3/4*b^2 + 3/2*b + 1/2
813            sage: (1+b) / (1-b) == c
814            True
815            sage: c * (1-b)
816            b + 1       
817        """
818        return self / other
819
820    def __nonzero__(self):
821        """
822        Return True if this number field element is nonzero.
823
824        EXAMPLES:
825            sage: m.<b> = CyclotomicField(17)
826            sage: m(0).__nonzero__()
827            False
828            sage: b.__nonzero__()
829            True
830
831        Nonzero is used by the bool command:
832            sage: bool(b + 1)
833            True
834            sage: bool(m(0))
835            False
836        """
837        return not IsZero_ZZX(self.__numerator)
838
839    cdef ModuleElement _neg_c_impl(self):
840        cdef NumberFieldElement x
841        x = self._new()
842        mul_ZZX_long(x.__numerator, self.__numerator, -1)
843        x.__denominator = self.__denominator
844        return x
845
846    def __int__(self):
847        """
848        Attempt to convert this number field element to a Python integer,
849        if possible.
850       
851        EXAMPLES:
852            sage: C.<I>=CyclotomicField(4)
853            sage: int(1/I)
854            Traceback (most recent call last):
855            ...
856            TypeError: cannot coerce nonconstant polynomial to int
857            sage: int(I*I)
858            -1
859
860            sage: K.<a> = NumberField(x^10 - x - 1)
861            sage: int(a)
862            Traceback (most recent call last):
863            ...
864            TypeError: cannot coerce nonconstant polynomial to int
865            sage: int(K(9390283))
866            9390283
867
868        The semantics are like in Python, so the value does not have
869        to preserved.
870            sage: int(K(393/29))
871            13           
872        """
873        return int(self.polynomial())
874
875    def __long__(self):
876        """
877        Attempt to convert this number field element to a Python long,
878        if possible.
879       
880        EXAMPLES:
881            sage: K.<a> = NumberField(x^10 - x - 1)
882            sage: long(a)
883            Traceback (most recent call last):
884            ...
885            TypeError: cannot coerce nonconstant polynomial to long
886            sage: long(K(1234))
887            1234L
888
889        The value does not have to be preserved, in the case of fractions.
890            sage: long(K(393/29))
891            13L       
892        """
893        return long(self.polynomial())
894
895    cdef void _parent_poly_c_(self, ZZX_c *num, ZZ_c *den):
896        raise NotImplementedError, "NumberFieldElement subclasses must override _parent_poly_c_()"
897        cdef long i
898        cdef ZZ_c coeff
899        cdef ntl_ZZX _num
900        cdef ntl_ZZ _den
901        if isinstance(self.parent(), sage.rings.number_field.number_field.NumberField_relative):
902            # ugly temp code
903            f = self.parent().absolute_polynomial()
904
905            __den = f.denominator()
906            (<Integer>ZZ(__den))._to_ZZ(den)
907
908            __num = f * __den
909            for i from 0 <= i <= __num.degree():
910                (<Integer>ZZ(__num[i]))._to_ZZ(&coeff)
911                SetCoeff( num[0], i, coeff )
912        else:
913            _num, _den = self.parent().polynomial_ntl()
914            num[0] = _num.x[0]
915            den[0] = _den.x[0]
916
917    cdef void _invert_c_(self, ZZX_c *num, ZZ_c *den):
918        """
919        Computes the numerator and denominator of the multiplicative inverse of this element.
920
921        Suppose that this element is x/d and the parent mod'ding polynomial is M/D.  The NTL function
922        XGCD( r, s, t, a, b ) computes r,s,t such that $r=s*a+t*b$.  We compute
923        XGCD( r, s, t, x*D, M*d ) and set
924        num=s*D*d
925        den=r
926
927        EXAMPLES:
928            I'd love to, but since we are dealing with c-types, I can't at this level.
929            Check __invert__ for doc-tests that rely on this functionality.
930        """
931        cdef ZZ_c parent_den
932        cdef ZZX_c parent_num
933        self._parent_poly_c_( &parent_num, &parent_den )
934
935        cdef ZZX_c t # unneeded except to be there
936        cdef ZZX_c a, b
937        mul_ZZX_ZZ( a, self.__numerator, parent_den )
938        mul_ZZX_ZZ( b, parent_num, self.__denominator )
939        XGCD_ZZX( den[0], num[0],  t, a, b, 1 )
940        mul_ZZX_ZZ( num[0], num[0], parent_den )
941        mul_ZZX_ZZ( num[0], num[0], self.__denominator )
942
943    def __invert__(self):
944        """
945        Returns the multiplicative inverse of self in the number field.
946
947        EXAMPLES:
948            sage: C.<I>=CyclotomicField(4)
949            sage: ~I
950            -I
951            sage: (2*I).__invert__()
952            -1/2*I
953        """
954        if IsZero_ZZX(self.__numerator):
955            raise ZeroDivisionError
956        cdef NumberFieldElement x
957        x = self._new()
958        self._invert_c_(&x.__numerator, &x.__denominator)
959        x._reduce_c_()
960        return x
961#        K = self.parent()
962#        quotient = K(1)._pari_('x') / self._pari_('x')
963#        if isinstance(K, sage.rings.number_field.number_field.NumberField_relative):
964#            return K(K.pari_rnf().rnfeltreltoabs(quotient))
965#        else:
966#            return K(quotient)
967
968    def _integer_(self):
969        """
970        Returns a rational integer if this element is actually a rational integer.
971
972        EXAMPLES:
973            sage: C.<I>=CyclotomicField(4)
974            sage: (~I)._integer_()
975            Traceback (most recent call last):
976            ...
977            TypeError: Unable to coerce -I to an integer
978            sage: (2*I*I)._integer_()
979            -2
980        """
981        if deg(self.__numerator) >= 1:
982            raise TypeError, "Unable to coerce %s to an integer"%self
983        return ZZ(self._rational_())
984
985    def _rational_(self):
986        """
987        Returns a rational number if this element is actually a rational number.
988
989        EXAMPLES:
990            sage: C.<I>=CyclotomicField(4)
991            sage: (~I)._rational_()
992            Traceback (most recent call last):
993            ...
994            TypeError: Unable to coerce -I to a rational
995            sage: (I*I/2)._rational_()
996            -1/2
997        """
998        if deg(self.__numerator) >= 1:
999            raise TypeError, "Unable to coerce %s to a rational"%self
1000        cdef Integer num
1001        num = PY_NEW(Integer)
1002        ZZX_getitem_as_mpz(&num.value, &self.__numerator, 0)
1003        return num / (<IntegerRing_class>ZZ)._coerce_ZZ(&self.__denominator)
1004
1005    def conjugate(self):
1006        """
1007        Return the complex conjugate of the number field element.  Currently,
1008        this is implemented for cyclotomic fields and quadratic extensions of Q. 
1009        It seems likely that there are other number fields for which the idea of
1010        a conjugate would be easy to compute.
1011
1012        EXAMPLES:
1013            sage: k.<I> = QuadraticField(-1)
1014            sage: I.conjugate()
1015            -I
1016            sage: (I/(1+I)).conjugate()
1017            -1/2*I + 1/2
1018            sage: z6=CyclotomicField(6).gen(0)
1019            sage: (2*z6).conjugate()
1020            -2*zeta6 + 2
1021
1022            sage: K.<b> = NumberField(x^3 - 2)
1023            sage: b.conjugate()
1024            Traceback (most recent call last):
1025            ...
1026            NotImplementedError: complex conjugation is not implemented (or doesn't make sense).
1027        """
1028        coeffs = self.parent().polynomial().list()
1029        if len(coeffs) == 3 and coeffs[2] == 1 and coeffs[1] == 0:
1030            # polynomial looks like x^2+d
1031            # i.e. we live in a quadratic extension of QQ
1032            if coeffs[0] > 0:
1033                gen = self.parent().gen()
1034                return self.polynomial()(-gen)
1035            else:
1036                return self
1037        elif isinstance(self.parent(), sage.rings.number_field.number_field.NumberField_cyclotomic):
1038            # We are in a cyclotomic field
1039            # Replace the generator zeta_n with (zeta_n)^(n-1)
1040            gen = self.parent().gen()
1041            return self.polynomial()(gen ** (gen.multiplicative_order()-1))
1042        else:
1043            raise NotImplementedError, "complex conjugation is not implemented (or doesn't make sense)."
1044
1045    def polynomial(self, var='x'):
1046        """
1047        Return the underlying polynomial corresponding to this
1048        number field element.
1049
1050        The resulting polynomial is currently *not* cached.
1051
1052        EXAMPLES:
1053            sage: K.<a> = NumberField(x^5 - x - 1)
1054            sage: f = (-2/3 + 1/3*a)^4; f
1055            1/81*a^4 - 8/81*a^3 + 8/27*a^2 - 32/81*a + 16/81
1056            sage: g = f.polynomial(); g
1057            1/81*x^4 - 8/81*x^3 + 8/27*x^2 - 32/81*x + 16/81
1058            sage: parent(g)
1059            Univariate Polynomial Ring in x over Rational Field
1060
1061        Note that the result of this function is not cached (should this
1062        be changed?):
1063            sage: g is f.polynomial()
1064            False       
1065        """
1066        return QQ[var](self._coefficients())
1067
1068    def _coefficients(self):
1069        """
1070        Return the coefficients of the underlying polynomial corresponding to this
1071        number field element.
1072
1073        OUTPUT:
1074             -- a list whose length corresponding to the degree of this element
1075                written in terms of a generator.
1076
1077        EXAMPLES:
1078           
1079        """
1080        coeffs = []
1081        cdef Integer den = (<IntegerRing_class>ZZ)._coerce_ZZ(&self.__denominator)
1082        cdef Integer numCoeff
1083        cdef int i
1084        for i from 0 <= i <= deg(self.__numerator):
1085            numCoeff = PY_NEW(Integer)
1086            ZZX_getitem_as_mpz(&numCoeff.value, &self.__numerator, i)
1087            coeffs.append( numCoeff / den )
1088        return coeffs
1089
1090    def denominator(self):
1091        """
1092        Return the denominator of this element, which is by definition
1093        the denominator of the corresponding polynomial
1094        representation.  I.e., elements of number fields are
1095        represented as a polynomial (in reduced form) modulo the
1096        modulus of the number field, and the denominator is the
1097        denominator of this polynomial.
1098
1099        EXAMPLES:
1100            sage: K.<z> = CyclotomicField(3)
1101            sage: a = 1/3 + (1/5)*z
1102            sage: print a.denominator()
1103            15
1104        """
1105        return (<IntegerRing_class>ZZ)._coerce_ZZ(&self.__denominator)
1106
1107    def _set_multiplicative_order(self, n):
1108        """
1109        Set the multiplicative order of this number field element.
1110
1111        WARNING -- use with caution -- only for internal use!  End
1112        users should never call this unless they have a very good
1113        reason to do so.
1114
1115        EXAMPLES:
1116            sage: K.<a> = NumberField(x^2 + x + 1)
1117            sage: a._set_multiplicative_order(3)
1118            sage: a.multiplicative_order()
1119            3
1120
1121        You can be evil with this so be careful.  That's why the function
1122        name begins with an underscore.
1123            sage: a._set_multiplicative_order(389)
1124            sage: a.multiplicative_order()
1125            389       
1126        """
1127        self.__multiplicative_order = n
1128
1129    def multiplicative_order(self):
1130        """
1131        Return the multiplicative order of this number field element.
1132
1133        EXAMPLES:
1134            sage: K.<z> = CyclotomicField(5)
1135            sage: z.multiplicative_order()
1136            5
1137            sage: (-z).multiplicative_order()
1138            10
1139            sage: (1+z).multiplicative_order()
1140            +Infinity       
1141        """
1142        if self.__multiplicative_order is not None:
1143            return self.__multiplicative_order
1144
1145        if self.is_rational_c():
1146            self.__multiplicative_order = self._rational_().multiplicative_order()
1147            return self.__multiplicative_order
1148
1149        if isinstance(self.parent(), sage.rings.number_field.number_field.NumberField_cyclotomic):
1150            t = self.parent()._multiplicative_order_table()
1151            f = self.polynomial()
1152            if t.has_key(f):
1153                self.__multiplicative_order = t[f]
1154                return self.__multiplicative_order
1155
1156        ####################################################################
1157        # VERY DUMB Algorithm to compute the multiplicative_order of
1158        # an element x of a number field K.
1159        #
1160        # 1. Find an integer B such that if n>=B then phi(n) > deg(K).
1161        #    For this use that for n>6 we have phi(n) >= log_2(n)
1162        #    (to see this think about the worst prime factorization
1163        #    in the multiplicative formula for phi.)
1164        # 2. Compute x, x^2, ..., x^B in order to determine the multiplicative_order.
1165        #
1166        # todo-- Alternative: Only do the above if we don't require an optional
1167        # argument which gives a multiple of the order, which is usually
1168        # something available in any actual application.
1169        #
1170        # BETTER TODO: Factor cyclotomic polynomials over K to determine
1171        # possible orders of elements?  Is there something even better?
1172        #
1173        ####################################################################
1174        d = self.parent().degree()
1175        B = max(7, 2**d+1)
1176        x = self
1177        i = 1
1178        while i < B:
1179            if x == 1:
1180                self.__multiplicative_order = i
1181                return self.__multiplicative_order
1182            x *= self
1183            i += 1
1184
1185        # it must have infinite order
1186        self.__multiplicative_order = sage.rings.infinity.infinity
1187        return self.__multiplicative_order
1188
1189    cdef bint is_rational_c(self):
1190        return deg(self.__numerator) == 0
1191               
1192    def trace(self):
1193        """
1194        Return the trace of this number field element.
1195
1196        EXAMPLES:
1197            sage: K.<a> = NumberField(x^3 -132/7*x^2 + x + 1); K
1198            Number Field in a with defining polynomial x^3 - 132/7*x^2 + x + 1
1199            sage: a.trace()
1200            132/7
1201            sage: (a+1).trace() == a.trace() + 3
1202            True       
1203        """
1204        K = self.parent().base_ring()
1205        return K(self._pari_('x').trace())
1206
1207    def norm(self):
1208        """
1209        Return the norm of this number field element.
1210
1211        EXAMPLES:
1212            sage: K.<a> = NumberField(x^3 + x^2 + x + -132/7); K
1213            Number Field in a with defining polynomial x^3 + x^2 + x - 132/7
1214            sage: a.norm()
1215            132/7
1216            sage: K(0).norm()
1217            0       
1218        """
1219        K = self.parent().base_ring()
1220        return K(self._pari_('x').norm())
1221
1222    def charpoly(self, var='x'):
1223        raise NotImplementedError, "Subclasses of NumberFieldElement must override charpoly()"
1224
1225    def minpoly(self, var='x'):
1226        """
1227        Return the minimal polynomial of this number field element.
1228
1229        EXAMPLES:
1230            sage: K.<a> = NumberField(x^2+3)
1231            sage: a.minpoly('x')
1232            x^2 + 3
1233            sage: R.<X> = K['X']
1234            sage: L.<b> = K.extension(X^2-(22 + a))
1235            sage: b.minpoly('t')
1236            t^4 + (-44)*t^2 + 487
1237            sage: b^2 - (22+a)
1238            0       
1239        """
1240        return self.charpoly(var).radical() # square free part of charpoly
1241
1242    def is_integral(self):
1243        r"""
1244        Determine if a number is in the ring of integers
1245        of this number field.
1246       
1247        EXAMPLES:
1248            sage: K.<a> = NumberField(x^2 + 23, 'a')
1249            sage: a.is_integral()
1250            True
1251            sage: t = (1+a)/2
1252            sage: t.is_integral()
1253            True
1254            sage: t.minpoly()
1255            x^2 - x + 6
1256            sage: t = a/2
1257            sage: t.is_integral()
1258            False
1259            sage: t.minpoly()
1260            x^2 + 23/4
1261        """
1262        return all([a in ZZ for a in self.minpoly()])
1263
1264    def matrix(self):
1265        r"""
1266        The matrix of right multiplication by the element on the power
1267        basis $1, x, x^2, \ldots, x^{d-1}$ for the number field.  Thus
1268        the {\em rows} of this matrix give the images of each of the $x^i$.
1269
1270        EXAMPLES:
1271
1272        Regular number field:
1273            sage: K.<a> = NumberField(QQ['x'].0^3 - 5)
1274            sage: M = a.matrix(); M
1275            [0 1 0]
1276            [0 0 1]
1277            [5 0 0]
1278            sage: M.base_ring() is QQ
1279            True
1280
1281        """
1282##         Relative number field:
1283##             sage: L.<b> = K.extension(K['x'].0^2 - 2)
1284##             sage: 1*b, b*b, b**3, b**6
1285##             (b, b^2, b^3, 6*b^4 - 10*b^3 - 12*b^2 - 60*b - 17)
1286##             sage: L.pari_rnf().rnfeltabstorel(b._pari_())
1287##             x - y
1288##             sage: L.pari_rnf().rnfeltabstorel((b**2)._pari_())
1289##             2
1290##             sage: M = b.matrix(); M
1291##             [0 1]
1292##             [3 0]
1293##             sage: M.base_ring() is K
1294##             True
1295
1296#         Absolute number field:
1297#             sage: M = L.absolute_field()[0].gen().matrix(); M
1298#             [  0   1   0   0   0   0]
1299#             [  0   0   1   0   0   0]
1300#             [  0   0   0   1   0   0]
1301#             [  0   0   0   0   1   0]
1302#             [  0   0   0   0   0   1]
1303#             [  2 -90 -27 -10   9   0]
1304#             sage: M.base_ring() is QQ
1305#             True
1306
1307#         More complicated relative number field:
1308#             sage: L.<b> = K.extension(K['x'].0^2 - a); L
1309#             Extension by x^2 + -a of the Number Field in a with defining polynomial x^3 - 5
1310#             sage: M = b.matrix(); M
1311#             [0 1]
1312#             [a 0]
1313#             sage: M.base_ring()
1314#             sage: M.base_ring() is K
1315#             True
1316        # Mutiply each power of field generator on
1317        # the left by this element; make matrix
1318        # whose rows are the coefficients of the result,
1319        # and transpose.
1320        if self.__matrix is None:
1321            K = self.parent()
1322            v = []
1323            x = K.gen()
1324            a = K(1)
1325            d = K.degree()
1326            for n in range(d):
1327                v += (a*self).list()
1328                a *= x
1329            k = K.base_ring()
1330            import sage.matrix.matrix_space
1331            M = sage.matrix.matrix_space.MatrixSpace(k, d)
1332            self.__matrix = M(v)
1333        return self.__matrix
1334
1335    def list(self):
1336        """
1337        Return list of coefficients of self written in terms of a power basis.
1338        """
1339        # Power basis list is total nonsense, unless the parent of self is an
1340        # absolute extension.
1341        raise NotImplementedError
1342       
1343       
1344cdef class NumberFieldElement_absolute(NumberFieldElement):
1345
1346    def _pari_(self, var='x'):
1347        """
1348        Return PARI C-library object corresponding to self.
1349
1350        EXAMPLES:
1351            sage: k.<j> = QuadraticField(-1)
1352            sage: j._pari_('j')
1353            Mod(j, j^2 + 1)
1354            sage: pari(j)
1355            Mod(x, x^2 + 1)
1356
1357            sage: y = QQ['y'].gen()
1358            sage: k.<j> = NumberField(y^3 - 2)
1359            sage: pari(j)
1360            Mod(x, x^3 - 2)
1361
1362        By default the variable name is 'x', since in PARI many variable
1363        names are reserved:
1364            sage: theta = polygen(QQ, 'theta')
1365            sage: M.<theta> = NumberField(theta^2 + 1)
1366            sage: pari(theta)
1367            Mod(x, x^2 + 1)
1368
1369        If you try do coerce a generator called I to PARI, hell may
1370        break loose:
1371            sage: k.<I> = QuadraticField(-1)
1372            sage: I._pari_('I')
1373            Traceback (most recent call last):
1374            ...
1375            PariError: forbidden (45)
1376
1377        Instead, request the variable be named different for the coercion:
1378            sage: pari(I)
1379            Mod(x, x^2 + 1)
1380            sage: I._pari_('i')
1381            Mod(i, i^2 + 1)
1382            sage: I._pari_('II')
1383            Mod(II, II^2 + 1)
1384        """
1385        try:
1386            return self.__pari[var]
1387        except KeyError:
1388            pass
1389        except TypeError:
1390            self.__pari = {}
1391        if var is None:
1392            var = self.parent().variable_name()
1393        f = self.polynomial()._pari_()
1394        gp = self.parent().polynomial()
1395        if gp.name() != 'x':
1396            gp = gp.change_variable_name('x')
1397        g = gp._pari_()
1398        gv = str(gp.parent().gen())
1399        if var != 'x':
1400            f = f.subst("x",var)
1401        if var != gv:
1402            g = g.subst(gv, var)
1403        h = f.Mod(g)
1404        self.__pari[var] = h
1405        return h
1406       
1407    cdef void _parent_poly_c_(self, ZZX_c *num, ZZ_c *den):
1408        cdef ntl_ZZX _num
1409        cdef ntl_ZZ _den
1410        _num, _den = self.parent().polynomial_ntl()
1411        num[0] = _num.x[0]
1412        den[0] = _den.x[0]
1413
1414    def charpoly(self, var='x'):
1415        r"""
1416        The characteristic polynomial of this element over $\Q$.
1417
1418        EXAMPLES:
1419
1420        We compute the charpoly of cube root of $2$.
1421
1422            sage: R.<x> = QQ[]
1423            sage: K.<a> = NumberField(x^3-2)
1424            sage: a.charpoly('x')
1425            x^3 - 2
1426
1427        """
1428        R = self.parent().base_ring()[var]
1429        return R(self._pari_('x').charpoly())
1430       
1431    def list(self):
1432        """
1433        Return list of coefficients of self written in terms of a power basis.
1434       
1435        EXAMPLE:
1436            sage: K.<z> = CyclotomicField(3)
1437            sage: (2+3/5*z).list()
1438            [2, 3/5]
1439            sage: (5*z).list()
1440            [0, 5]
1441            sage: K(3).list()
1442            [3, 0]
1443        """
1444        n = self.parent().degree()
1445        v = self._coefficients()
1446        z = sage.rings.rational.Rational(0)
1447        return v + [z]*(n - len(v))
1448       
1449       
1450cdef class NumberFieldElement_quadratic(NumberFieldElement_absolute):
1451    # (a + b sqrt(disc)) / denom
1452    def __init__(self, parent, f):
1453        NumberFieldElement_absolute.__init__(self, parent, f) # this is heavy, can I avoid it for f in Rational, tuple?
1454        self.disc = parent.discriminant()
1455        cdef Integer a, b, denom
1456        cdef Rational ad, bd
1457       
1458        cdef Py_ssize_t fdeg
1459        cdef mpz_t tmp
1460        cdef NumberFieldElement_quadratic gen
1461       
1462        if PY_TYPE_CHECK_EXACT(f, Rational):
1463            mpz_set(self.a, mpq_numref((<Rational>f).value))
1464            mpz_set_ui(self.b, 0)
1465            mpz_set(self.denom, mpq_denref((<Rational>f).value))
1466           
1467        elif PY_TYPE_CHECK_EXACT(f, tuple) and len(f) == 2:
1468            ad, bd = f
1469            mpz_set(self.a, a.value)
1470            mpz_set(self.b, b.value)
1471            mpz_lcm(self.denom, mpq_denref(ad.value), mpq_denref(bd.value))
1472            mpz_divexact(self.a, self.denom, mpq_denref(ad.value))
1473            mpz_mul(self.a, self.a, mpq_numref(ad.value))
1474            mpz_divexact(self.b, self.denom, mpq_denref(bd.value))
1475            mpz_mul(self.b, self.b, mpq_numref(bd.value))
1476
1477        elif PY_TYPE_CHECK_EXACT(f, tuple) and len(f) == 3:
1478            a, b, denom = f
1479            mpz_set(self.a, a.value)
1480            mpz_set(self.b, b.value)
1481            mpz_set(self.denom, denom.value)
1482            self._reduce_c()
1483           
1484        else:
1485            fdeg = deg(self.__numerator) # poly in generator
1486            if fdeg > -1:
1487                ZZX_getitem_as_mpz(&self.a,  &self.__numerator, 0)
1488                if fdeg > 0:
1489                    gen = parent.gen() # should this be cached?
1490                    mpz_init(tmp)
1491                    ZZX_getitem_as_mpz(&tmp,  &self.__numerator, 1)
1492                    mpz_mul(self.b, tmp, gen.a) # using self.b as temp
1493                    mpz_add(self.a, self.a, self.b)
1494                    mpz_mul(self.b, tmp, gen.b) # really setting self.b
1495                    mpz_clear(tmp)
1496                else:
1497                    mpz_set_ui(self.b, 0)
1498                denom = (<IntegerRing_class>ZZ)._coerce_ZZ(&self.__denominator)
1499                mpz_set(self.denom, denom.value)
1500            else:
1501                mpz_set_ui(self.a, 0)
1502                mpz_set_ui(self.b, 0)
1503                mpz_set_ui(self.denom, 1)
1504
1505    def __new__(self):
1506        # we do NOT own self.disc.value
1507        mpz_init(self.a)
1508        mpz_init(self.b)
1509        mpz_init(self.denom)
1510       
1511    def __dealloc__(self):
1512        # we do NOT own self.disc.value
1513        mpz_clear(self.a)
1514        mpz_clear(self.b)
1515        mpz_clear(self.denom)
1516   
1517    def parts(self):
1518        cdef Rational ad = <Rational>PY_NEW(Rational), bd = <Rational>PY_NEW(Rational)
1519        mpz_set(mpq_numref(ad.value), self.a)
1520        mpz_set(mpq_denref(ad.value), self.denom)
1521        mpz_set(mpq_numref(bd.value), self.b)
1522        mpz_set(mpq_denref(bd.value), self.denom)
1523        return ad, bd
1524       
1525    cdef bint is_sqrt_disc(self):
1526        return mpz_cmp_ui(self.denom, 1)==0 and mpz_cmp_ui(self.a,0)==0 and mpz_cmp_ui(self.b, 1)==0
1527               
1528               
1529#########################################################
1530# Arithmatic
1531#########################################################
1532
1533    cdef void _reduce_c_(self):
1534        cdef mpz_t gcd
1535        mpz_init(gcd)
1536        mpz_gcd(gcd, self.a, self.b)
1537        mpz_gcd(gcd, gcd, self.denom)
1538        if mpz_cmp_si(gcd, 1): # != 0 (i.e. it is not 1)
1539            mpz_divexact(self.a, self.a, gcd)
1540            mpz_divexact(self.b, self.b, gcd)
1541            mpz_divexact(self.denom, self.denom, gcd)
1542        mpz_clear(gcd)
1543       
1544    cdef ModuleElement _add_c_impl(self, ModuleElement other_m):
1545        cdef NumberFieldElement_quadratic other = <NumberFieldElement_quadratic>other_m
1546        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1547        res.disc = self.disc
1548        mpz_add(res.a, self.a, other.a)
1549        mpz_add(res.b, self.b, other.b)
1550        mpz_set(res.denom, self.denom)
1551        res._reduce_c_()
1552        return res
1553   
1554    cdef ModuleElement _sub_c_impl(self, ModuleElement other_m):
1555        cdef NumberFieldElement_quadratic other = <NumberFieldElement_quadratic>other_m
1556        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1557        res.disc = self.disc
1558        mpz_sub(res.a, self.a, other.a)
1559        mpz_sub(res.b, self.b, other.b)
1560        mpz_set(res.denom, self.denom)
1561        res._reduce_c_()
1562        return res
1563       
1564    def __neg__(self):
1565        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1566        res.disc = self.disc
1567        mpz_neg(res.a, self.a)
1568        mpz_neg(res.b, self.b)
1569        mpz_set(res.denom, self.denom)
1570        return res
1571   
1572    cdef RingElement _mul_c_impl(self, RingElement other_m):
1573        cdef NumberFieldElement_quadratic other = <NumberFieldElement_quadratic>other_m
1574        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1575        res.disc = self.disc
1576        cdef mpz_t tmp
1577        mpz_init(tmp)
1578       
1579        if mpz_size(self.a) < 2: # could I use a macro instead?
1580            # Do it the traditional way
1581            mpz_mul(res.a, self.a, other.a)
1582            mpz_mul(tmp, self.b, other.b)
1583            mpz_mul(tmp, tmp, self.disc.value)
1584            mpz_add(res.a, res.a, tmp)
1585           
1586            mpz_mul(res.b, self.a, other.b)
1587            mpz_mul(tmp, self.a, other.b)
1588            mpz_add(res.b, res.b, tmp)
1589           
1590        else:
1591            # Karatsuba
1592            mpz_add(res.a, self.a, other.a) # using res.a as tmp
1593            mpz_add(tmp, self.b, other.b)
1594            mpz_mul(res.b, res.a, tmp) # res.b = (self.a + other.a)(self.b + other.b)
1595           
1596            mpz_mul(res.a, self.a, other.a)
1597            mpz_sub(res.b, res.b, res.a)
1598            mpz_mul(tmp, self.b, other.b)
1599            mpz_sub(res.b, res.b, tmp)
1600            mpz_mul(tmp, tmp, self.disc.value)
1601            mpz_add(res.a, tmp, tmp)
1602       
1603        mpz_clear(tmp)
1604       
1605        mpz_mul(res.denom, self.denom, other.denom)
1606        res._reduce_c_()
1607
1608        return res
1609       
1610    cdef ModuleElement _rmul_c_impl(self, RingElement _c):
1611        cdef Rational c = <Rational>_c
1612        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1613        res.disc = self.disc
1614        mpz_mul(res.a, self.a, mpq_numref(c.value))
1615        mpz_mul(res.b, self.b, mpq_numref(c.value))
1616        mpz_mul(res.denom, self.denom, mpq_denref(c.value))
1617        res._reduce_c_()
1618        return res
1619       
1620    cdef ModuleElement _lmul_c_impl(self, RingElement _c):
1621        cdef Rational c = <Rational>_c
1622        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1623        res.disc = self.disc
1624        mpz_mul(res.a, self.a, mpq_numref(c.value))
1625        mpz_mul(res.b, self.b, mpq_numref(c.value))
1626        mpz_mul(res.denom, self.denom, mpq_denref(c.value))
1627        res._reduce_c_()
1628        return res
1629       
1630    cdef RingElement _div_c_impl(self, RingElement other):
1631        return self * ~other
1632       
1633    def __invert__(self):
1634        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1635        res.disc = self.disc
1636        cdef mpz_t tmp, gcd
1637        mpz_init(tmp)
1638        mpz_init(gcd)
1639       
1640        mpz_gcd(gcd, self.a, self.b)
1641        if mpz_cmp_si(gcd, 1): # != 0 (i.e. it is not 1)
1642            # cancel out g (g(a'-b'd)) / (g^2 (a'^2-b'^2d^2))
1643            mpz_divexact(res.a, self.a, gcd)
1644            mpz_divexact(res.b, self.b, gcd)
1645            mpz_neg(res.b, res.b)
1646        else:
1647            mpz_set(res.a, self.a)
1648            mpz_neg(res.b, self.b)
1649
1650        mpz_pow_ui(res.denom, res.a, 2)
1651        mpz_pow_ui(tmp, res.b, 2)
1652        mpz_mul(tmp, tmp, self.disc.value)
1653        mpz_sub(res.denom, res.denom, tmp)
1654        # need to multiply the leftover g back in
1655        mpz_mul(res.denom, res.denom, gcd)
1656           
1657        mpz_mul(res.denom, res.denom, self.denom)
1658
1659        mpz_clear(tmp)
1660        mpz_clear(gcd)
1661
1662        res._reduce_c_()
1663        return res
1664       
1665    cdef NumberFieldElement conjugate_c(self):
1666        cdef NumberFieldElement_quadratic res = <NumberFieldElement_quadratic>self._new_c()
1667        res.disc = self.disc
1668        mpz_set(res.a, self.a)
1669        mpz_neg(res.b, self.b)
1670        mpz_set(res.denom, self.denom)
1671        return res
1672       
1673#################################################################################
1674# We must override everything that makes uses of self.__numerator/__denominator
1675#################################################################################
1676           
1677    cdef int _cmp_c_impl(self, sage.structure.element.Element _right) except -2:
1678        cdef NumberFieldElement_quadratic right = _right
1679        return not mpz_cmp(self.a, right.a)==0  \
1680            or not mpz_cmp(self.b, right.b)==0  \
1681            or not mpz_cmp(self.denom, right.denom) == 0
1682
1683    def __nonzero__(self):
1684        return mpz_cmp_ui(self.a, 0) != 0 or mpz_cmp_ui(self.b, 0) != 0
1685   
1686    def _integer_(self):
1687        cdef Integer res
1688        if mpz_cmp_ui(self.b, 0) != 0 or mpz_cmp_ui(self.denom, 1) != 0:
1689            raise TypeError, "Unable to coerce %s to an integer"%self
1690        else:
1691            res = PY_NEW(Integer)
1692            mpz_set(res.value, self.a)
1693            return res
1694
1695    def _rational_(self):
1696        cdef Rational res
1697        if mpz_cmp_ui(self.b, 0)!=0:
1698            raise TypeError, "Unable to coerce %s to a rational"%self
1699        else:
1700            res = <Rational>PY_NEW(Rational)
1701            mpz_set(mpq_numref(res.value), self.a)
1702            mpz_set(mpq_denref(res.value), self.denom)
1703            mpq_canonicalize(res.value)
1704            return res
1705
1706    def _coefficients(self):
1707        # In terms of the generator...
1708        cdef NumberFieldElement_quadratic gen = self._parent.gen() # should this be cached?
1709        cdef Rational const = <Rational>PY_NEW(Rational), lin = <Rational>PY_NEW(Rational)
1710        if gen.is_sqrt_disc():
1711            # gen = sqrt(disc)
1712            ad, bd = self.parts()
1713            return [ad,bd]
1714        else:
1715            alpha, beta = gen.parts()
1716            scale = bd/beta
1717            return [ad - scale*alpha, scale]
1718       
1719    def denominator(self):
1720        cdef Integer denom = PY_NEW(Integer)
1721        mpz_set(denom.value, self.denom)
1722        return denom
1723       
1724    cdef bint is_rational_c(self):
1725        return mpz_cmp_ui(self.b, 0) == 0
1726
1727#########################################################
1728# Some things are so much easier to compute
1729#########################################################
1730
1731    def trace(self):
1732        # trace = 2*a
1733        cdef Rational res = <Rational>PY_NEW(Rational)
1734        if mpz_odd_p(self.denom):
1735            mpz_mul_2exp(mpq_numref(res.value), self.a, 1)
1736            mpz_set(mpq_denref(res.value), self.denom)
1737        else:
1738            mpz_set(mpq_numref(res.value), self.a)
1739            mpz_divexact_ui(mpq_denref(res.value), self.denom, 2)
1740
1741    def norm(self):
1742        # norm = a^2 - d b^2
1743        cdef Rational res = <Rational>PY_NEW(Rational)
1744        mpz_pow_ui(mpq_numref(res.value), self.a, 2)
1745        mpz_pow_ui(mpq_denref(res.value), self.b, 2) # use as temp
1746        mpz_mul(mpq_denref(res.value), mpq_denref(res.value), self.disc.value)
1747        mpz_sub(mpq_numref(res.value), mpq_numref(res.value), mpq_denref(res.value))
1748        mpz_pow_ui(mpq_denref(res.value), self.denom, 2)
1749        mpq_canonicalize(res.value)
1750        return res
1751       
1752    def charpoly(self, var='x'):
1753        r"""
1754        The characteristic polynomial of this element over $\Q$.
1755
1756        EXAMPLES:
1757
1758        We compute the charpoly of cube root of $2$.
1759
1760            sage: R.<x> = QQ[]
1761            sage: K.<a> = NumberField(x^2-7)
1762            sage: a.charpoly('x')
1763            x^2 - 7
1764
1765        """
1766        R = QQ[var]
1767        return R([self.norm(), -self.trace(), 1])
1768       
1769
1770cdef class NumberFieldElement_relative(NumberFieldElement):
1771
1772    def _pari_(self, var='x'):
1773        """
1774        Return PARI C-library object corresponding to self.
1775
1776        EXAMPLES:
1777        By default the variable name is 'x', since in PARI many variable
1778        names are reserved.
1779            sage: y = QQ['y'].gen()
1780            sage: k.<j> = NumberField([y^3 - 2, y^2 - 7])
1781            sage: pari(j)
1782            Mod(42/5515*x^5 - 9/11030*x^4 - 196/1103*x^3 + 273/5515*x^2 + 10281/5515*x + 4459/11030, x^6 - 21*x^4 + 4*x^3 + 147*x^2 + 84*x - 339)
1783            sage: j^2
1784            7
1785            sage: pari(j)^2
1786            Mod(7, x^6 - 21*x^4 + 4*x^3 + 147*x^2 + 84*x - 339)
1787        """
1788        try:
1789            return self.__pari[var]
1790        except KeyError:
1791            pass
1792        except TypeError:
1793            self.__pari = {}
1794        if var is None:
1795            var = self.parent().variable_name()
1796        f = self.polynomial()._pari_()
1797        g = str(self.parent().pari_polynomial())
1798        base = self.parent().base_ring()
1799        gsub = base.gen()._pari_()
1800        gsub = str(gsub).replace('x', 'y')
1801        g = g.replace('y', gsub)
1802        h = f.Mod(g)
1803        self.__pari[var] = h
1804        return h
1805       
1806    cdef void _parent_poly_c_(self, ZZX_c *num, ZZ_c *den):
1807        cdef long i
1808        cdef ZZ_c coeff
1809        cdef ntl_ZZX _num
1810        cdef ntl_ZZ _den
1811        # ugly temp code
1812        f = self.parent().absolute_polynomial()
1813
1814        __den = f.denominator()
1815        (<Integer>ZZ(__den))._to_ZZ(den)
1816
1817        __num = f * __den
1818        for i from 0 <= i <= __num.degree():
1819            (<Integer>ZZ(__num[i]))._to_ZZ(&coeff)
1820            SetCoeff( num[0], i, coeff )
1821
1822    def __repr__(self):
1823        K = self.parent()
1824        # Compute representation of self in terms of relative vector space.
1825        w = self.vector()
1826        R = K.base_field()[K.variable_name()]
1827        return repr(R(w.list()))
1828
1829    def vector(self):
1830        return self.parent().vector_space()[2](self)
1831
1832    def charpoly(self, var='x'):
1833        r"""
1834        The characteristic polynomial of this element over $\Q$.
1835
1836        EXAMPLES:
1837
1838        We construct a relative extension and find the characteristic
1839        polynomial over $\Q$.
1840
1841            sage: R.<x> = QQ[]
1842            sage: K.<a> = NumberField(x^3-2)
1843            sage: S.<X> = K[]
1844            sage: L.<b> = NumberField(X^3 + 17); L
1845            Number Field in b with defining polynomial X^3 + 17 over its base field
1846            sage: b.charpoly ()
1847            x^9 + 51*x^6 + 867*x^3 + 4913
1848            sage: b.charpoly()(b)
1849            0
1850            sage: a = L.0; a
1851            b
1852            sage: a.charpoly('x')
1853            x^9 + 51*x^6 + 867*x^3 + 4913
1854            sage: a.charpoly('y')
1855            y^9 + 51*y^6 + 867*y^3 + 4913
1856        """
1857        R = self.parent().base_ring()[var]
1858        g = self.polynomial()  # in QQ[x]
1859        f = self.parent().pari_polynomial()  # # field is QQ[x]/(f)
1860        return R( (g._pari_().Mod(f)).charpoly() )
1861
1862## This might be useful for computing relative charpoly.
1863## BUT -- currently I don't even know how to view elements
1864## as being in terms of the right thing, i.e., this code
1865## below as is lies.
1866##             nf = self.parent()._pari_base_nf()
1867##             prp = self.parent().pari_relative_polynomial()
1868##             elt = str(self.polynomial()._pari_())
1869##             return R(nf.rnfcharpoly(prp, elt))
1870##         # return self.matrix().charpoly('x')
1871
1872
1873cdef class OrderElement_absolute(NumberFieldElement_absolute):
1874    """
1875    Element of an order in an absolute number field.
1876
1877    EXAMPLES:
1878        sage: k.<a> = NumberField(x^2 + 1)
1879    """
1880    def __init__(self, order, f):
1881        K = order.number_field()
1882        NumberFieldElement_absolute.__init__(self, K, f)
1883        self._order = order
1884
1885cdef class OrderElement_relative(NumberFieldElement_relative):
1886    """
1887    Element of an order in a relative number field.
1888    """
1889    def __init__(self, order, f):
1890        K = order.number_field()
1891        NumberFieldElement_relative.__init__(self, K, f)
1892        self._order = order
1893
1894
1895
1896
1897class CoordinateFunction:
1898    def __init__(self, alpha, W, to_V):
1899        self.__alpha = alpha
1900        self.__W = W
1901        self.__to_V = to_V
1902        self.__K = alpha.parent()
1903
1904    def __repr__(self):
1905        return "Coordinate function that writes elements in terms of the powers of %s"%self.__alpha
1906
1907    def alpha(self):
1908        return self.__alpha
1909
1910    def __call__(self, x):
1911        return self.__W.coordinates(self.__to_V(self.__K(x)))
1912   
Note: See TracBrowser for help on using the repository browser.