Ticket #11890: 11890_rebased.patch

File 11890_rebased.patch, 13.6 KB (added by jdemeyer, 9 years ago)
  • sage/libs/pari/gen.pyx

    # HG changeset patch
    # User Jeroen Demeyer <jdemeyer@cage.ugent.be>
    # Date 1317743222 -7200
    # Node ID 3e42a268abfbe60dc0950ab2161042184122aa94
    # Parent  4b285df55c4c596dbea3c62a603a5c5261e9a3ff
    Add "important" argument to pari_nf() to allow factoring of polynomials over big number fields.
    
    diff --git a/sage/libs/pari/gen.pyx b/sage/libs/pari/gen.pyx
    a b  
    446446        return self.g
    447447   
    448448    def list(self):
    449         if typ(self.g) == t_POL:
    450             raise NotImplementedError, \
    451                 "please report, t_POL.list() is broken, should not be used!"
    452         if typ(self.g) == t_SER:
    453             raise NotImplementedError, \
    454                 "please report, t_SER.list() is broken, should not be used!"
    455         #return list(self.Vecrev())
     449        """
     450        Convert self to a list of PARI gens.
     451
     452        EXAMPLES:
     453
     454        A PARI vector becomes a Sage list::
     455
     456            sage: L = pari("vector(10,i,i^2)").list()
     457            sage: L
     458            [1, 4, 9, 16, 25, 36, 49, 64, 81, 100]
     459            sage: type(L)
     460            <type 'list'>
     461            sage: type(L[0])
     462            <type 'sage.libs.pari.gen.gen'>
     463
     464        For polynomials, list() behaves as for ordinary Sage polynomials::
     465
     466            sage: pol = pari("x^3 + 5/3*x"); pol.list()
     467            [0, 5/3, 0, 1]
     468
     469        For power series, we get all coefficients, including leading and
     470        trailing zeros::
     471
     472            sage: R.<x> = LaurentSeriesRing(QQ)
     473            sage: s = x^2 + O(x^8)
     474            sage: s.list()
     475            [1]
     476            sage: pari(s).list()
     477            [0, 0, 1, 0, 0, 0, 0, 0]
     478
     479        For matrices, we get a list of columns::
     480
     481            sage: M = matrix(ZZ,3,2,[1,4,2,5,3,6]); M
     482            [1 4]
     483            [2 5]
     484            [3 6]
     485            sage: pari(M).list()
     486            [[1, 2, 3]~, [4, 5, 6]~]
     487
     488        For "scalar" types, we get a 1-element list containing ``self``::
     489
     490            sage: pari("42").list()
     491            [42]
     492        """
     493        if typ(self.g) == t_POL or typ(self.g) == t_SER:
     494            return list(self.Vecrev())
    456495        return list(self.Vec())
    457496       
    458497    def __reduce__(self):
     
    75027541    def __call__(self, x):
    75037542        return self.eval(x)
    75047543
     7544    def factornf(self, t):
     7545        """
     7546        Factorization of the polynomial ``self`` over the number field
     7547        defined by the polynomial ``t``.  This does not require that `t`
     7548        is integral, nor that the discriminant of the number field can be
     7549        factored.
     7550
     7551        EXAMPLES::
     7552
     7553            sage: x = polygen(QQ)
     7554            sage: K.<a> = NumberField(x^2 - 1/8)
     7555            sage: pari(x^2 - 2).factornf(K.pari_polynomial("a"))
     7556            [x + Mod(-4*a, 8*a^2 - 1), 1; x + Mod(4*a, 8*a^2 - 1), 1]
     7557        """
     7558        t0GEN(t)
     7559        sig_on()
     7560        return self.new_gen(polfnf(self.g, t0))
     7561
    75057562    def factorpadic(self, p, long r=20, long flag=0):
    75067563        """
    75077564        self.factorpadic(p,r=20,flag=0): p-adic factorization of the
  • sage/rings/number_field/number_field.py

    diff --git a/sage/rings/number_field/number_field.py b/sage/rings/number_field/number_field.py
    a b  
    26142614            self.__pari_polynomial = polypari
    26152615            return self.__pari_polynomial
    26162616
    2617     def pari_nf(self):
     2617    def pari_nf(self, important=True):
    26182618        """
    26192619        PARI number field corresponding to this field.
    26202620       
    26212621        This is the number field constructed using nfinit(). This is the same
    26222622        as the number field got by doing pari(self) or gp(self).
     2623
     2624        INPUT:
     2625
     2626        - ``important`` -- (default: True) bool.  If False, raise a
     2627          ``RuntimeError`` if we need to do a difficult discriminant
     2628          factorization.  Useful when the PARI nf structure is useful
     2629          but not strictly required, such as for factoring polynomials
     2630          over this number field.
    26232631       
    26242632        EXAMPLES::
    26252633       
     
    26472655            ...
    26482656            TypeError: Unable to coerce number field defined by non-integral polynomial to PARI.
    26492657
    2650         We illustrate the maximize_at_primes and assume_disc_small parameters::
    2651        
    2652             sage: p = next_prime(10^40); q = next_prime(p)
    2653 
    2654         The following would take a very long time without the
    2655         maximize_at_primes option above::
     2658        With ``important=False``, we simply bail out if we cannot
     2659        easily factor the discriminant::
     2660
     2661            sage: p = next_prime(10^40); q = next_prime(10^41)
     2662            sage: K.<a> = NumberField(x^2 - p*q)
     2663            sage: K.pari_nf(important=False)
     2664            Traceback (most recent call last):
     2665            ...
     2666            RuntimeError: Unable to factor discriminant with trial division
     2667
     2668        Next, we illustrate the ``maximize_at_primes`` and ``assume_disc_small``
     2669        parameters of the NumberField contructor. The following would take a
     2670        very long time without the ``maximize_at_primes`` option::
    26562671       
    26572672            sage: K.<a> = NumberField(x^2 - p*q, maximize_at_primes=[p])
    26582673            sage: K.pari_nf()
     
    26702685            return self.__pari_nf
    26712686        except AttributeError:
    26722687            f = self.pari_polynomial()
    2673             self.__pari_nf = pari([f, self._pari_integral_basis()]).nfinit()
     2688            self.__pari_nf = pari([f, self._pari_integral_basis(important=important)]).nfinit()
    26742689            return self.__pari_nf
    26752690
    26762691    def pari_zk(self):
     
    40144029        """
    40154030        return self.maximal_order(v=v).basis()
    40164031   
    4017     def _pari_integral_basis(self, v=None):
     4032    def _pari_integral_basis(self, v=None, important=True):
    40184033        """
    40194034        Internal function returning an integral basis of this number field in
    40204035        PARI format.
    40214036       
    40224037        INPUT:
    40234038       
    4024        
    4025         -  ``v`` - None, a prime, or a list of primes. See the
     4039        -  ``v`` -- None, a prime, or a list of primes. See the
    40264040           documentation for self.maximal_order.
    40274041       
     4042        - ``important`` -- (default:True) bool.  If False, raise a
     4043          ``RuntimeError`` if we need to do a difficult discriminant
     4044          factorization.  Useful when the PARI nf structure is useful
     4045          but not strictly required.
    40284046       
    40294047        EXAMPLES::
    40304048       
     
    40504068        v = self._normalize_prime_list(v)
    40514069        try:
    40524070            return self._integral_basis_dict[v]
    4053         except:
     4071        except (AttributeError, KeyError):
    40544072            f = self.pari_polynomial()
    4055             if len(v) == 0:
    4056                 B = f.nfbasis(1 if self._assume_disc_small else 0)
    4057             else:
     4073            if len(v) > 0:
    40584074                m = self._pari_disc_factorization_matrix(v)
    40594075                B = f.nfbasis(fa = m)
     4076            elif self._assume_disc_small:
     4077                B = f.nfbasis(1)
     4078            elif not important:
     4079                # Trial divide the discriminant
     4080                m = self.pari_polynomial().poldisc().abs().factor(limit=0)
     4081                # Since we only need a *squarefree* factorization, we need
     4082                # trial division up to D^(1/3) instead of D^(1/2).
     4083                trialdivlimit = pari(pari._primelimit()**3)
     4084                if all([ p < trialdivlimit or p.isprime() for p in m[0] ]):
     4085                    B = f.nfbasis(fa = m)
     4086                else:
     4087                    raise RuntimeError, "Unable to factor discriminant with trial division"
     4088            else:
     4089                B = f.nfbasis()
    40604090
    40614091            self._integral_basis_dict[v] = B
    40624092            return B
     
    64536483              To:   Number Field in a0 with defining polynomial x^2 - 1/2*c0 over its base field
    64546484              Defn: z |--> a0)
    64556485
    6456         We can relativize over a really large field.  This requires great care
    6457         to not factor or do any operation that would trigger a pari nfinit()
    6458         internally::
     6486        We can relativize over a really large field::
    64596487
    64606488            sage: K.<a> = CyclotomicField(3^3*2^3)
    64616489            sage: R = K.relativize(a^(3^2), 't'); R
     
    65306558        # the generator of the intermediate field.
    65316559
    65326560        # this strange step avoids computing an embedding of the base_field L
    6533         # into self, a computation which triggers a pari nfinit().  Without
    6534         # this, the relativize done over a huge field above is not feasible.
     6561        # into self.
    65356562        base_hom = L.hom([alpha], self)
    65366563        from_M = M.Hom(self)([self.gen()], base_hom=base_hom, check=True)
    65376564
     
    75607587            if p % n == 1:
    75617588                return p
    75627589   
    7563     def _pari_integral_basis(self, v=None):
     7590    def _pari_integral_basis(self, v=None, important=True):
    75647591        """
    75657592        Internal function returning an integral basis of this number field in
    75667593        PARI format.
    75677594       
    75687595        This field is cyclomotic, so this is a trivial computation, since
    7569         the power basis on the generator is an integral basis. Thus the v
    7570         parameter is ignored.
     7596        the power basis on the generator is an integral basis. Thus the ``v``
     7597        and ``important`` parameters are ignored.
    75717598       
    75727599        EXAMPLES::
    75737600       
  • sage/rings/polynomial/polynomial_element.pyx

    diff --git a/sage/rings/polynomial/polynomial_element.pyx b/sage/rings/polynomial/polynomial_element.pyx
    a b  
    27232723            sage: pol.factor()
    27242724            (t - 2*a^3 + a) * (t - 4/3*a^3 + 2/3*a) * (t - 2/3*a^3 + 1/3*a)
    27252725
     2726        Test that this factorization really uses ``nffactor()`` internally::
     2727
     2728            sage: pari.default("debug", 3)
     2729            sage: F = pol.factor()
     2730
     2731            Entering nffactor:
     2732            ...
     2733            sage: pari.default("debug", 0)
     2734
    27262735        Test that ticket #10369 is fixed::
    27272736
    27282737            sage: x = polygen(QQ)
     
    27412750            sage: pol.factor()
    27422751            x^5 * (x^5 + (4/7*a - 6/7)*x^4 + (9/49*a^2 - 3/7*a + 15/49)*x^3 + (8/343*a^3 - 32/343*a^2 + 40/343*a - 20/343)*x^2 + (5/2401*a^4 - 20/2401*a^3 + 40/2401*a^2 - 5/343*a + 15/2401)*x - 6/16807*a^4 + 12/16807*a^3 - 18/16807*a^2 + 12/16807*a - 6/16807)
    27432752
     2753        Factoring over a number field over which we cannot factor the
     2754        discriminant::
     2755
     2756            sage: p = next_prime(10^50); q = next_prime(10^51); n = p*q;
     2757            sage: K.<a> = QuadraticField(p*q)
     2758            sage: R.<x> = PolynomialRing(K)
     2759            sage: factor(x^2 + 1)
     2760            x^2 + 1
     2761            sage: factor( (x - a) * (x + 2*a) )
     2762            (x - a) * (x + 2*a)
    27442763        """
    27452764        # PERFORMANCE NOTE:
    27462765        #     In many tests with SMALL degree PARI is substantially
     
    28402859        elif is_NumberField(R):
    28412860       
    28422861            if R.degree() == 1:
    2843                
    28442862                factors = self.change_ring(QQ).factor()
    28452863                return Factorization([(self._parent(p), e) for p, e in factors], R(factors.unit()))
    28462864       
    2847             if R.defining_polynomial().denominator() == 1:
    2848                 v = [ x._pari_("a") for x in self.list() ]
    2849                 f = pari(v).Polrev()
    2850                 Rpari = R.pari_nf()
    2851                 if (Rpari.variable() != "a"):
    2852                     Rpari = copy.copy(Rpari)
    2853                     Rpari[0] = Rpari[0]("a")
    2854                     Rpari[6] = [ x("a") for x in Rpari[6] ]
     2865            # Use "a" for the number field variable
     2866            v = [ c._pari_("a") for c in self.list() ]
     2867            f = pari(v).Polrev()
     2868            try:
     2869                # Try to compute the PARI nf structure with important=False.
     2870                # This will raise RuntimeError if the computation is too
     2871                # difficult.  It will raise TypeError if the defining
     2872                # polynomial is not integral.
     2873                Rpari = R.pari_nf(important=False).nf_subst("'a")
     2874                # Factor using nffactor()
    28552875                G = list(Rpari.nffactor(f))
    2856                 # PARI's nffactor() ignores the unit, _factor_pari_helper()
    2857                 # adds back the unit of the factorization.
    2858                 return self._factor_pari_helper(G)
    2859 
    2860             else:
    2861 
    2862                 Rdenom = R.defining_polynomial().denominator()
    2863 
    2864                 new_Rpoly = (R.defining_polynomial() * Rdenom).change_variable_name("a")
    2865 
    2866                 Rpari, Rdiff = new_Rpoly._pari_().nfinit(3)
    2867 
    2868                 AZ = QQ['z']
    2869                 Raux = NumberField(AZ(Rpari[0]),'alpha')
    2870 
    2871                 S, gSRaux, fRauxS = Raux.change_generator(Raux(Rdiff))
    2872 
    2873                 phi_RS = R.Hom(S)([S.gen(0)])
    2874                 phi_SR = S.Hom(R)([R.gen(0)])
    2875 
    2876                 unit = self.leading_coefficient()
    2877                 temp_f = self * 1/unit
    2878 
    2879                 v = [ gSRaux(phi_RS(x))._pari_("a") for x in temp_f.list() ]
    2880                 f = pari(v).Polrev()
    2881 
    2882                 pari_factors = Rpari.nffactor(f)
    2883 
    2884                 factors = [ ( self.parent([ phi_SR(fRauxS(Raux(pari_factors[0][i][j])))
    2885                                             for j in range(len(pari_factors[0][i])) ]) ,
    2886                              int(pari_factors[1][i]) )
    2887                             for i in range(pari_factors.nrows()) ]
    2888 
    2889                 return Factorization(factors, unit)
    2890        
     2876            except (RuntimeError, TypeError):
     2877                # Use factornf() which only needs the defining polynomial
     2878                # and which does not require an integral polynomial.
     2879                G = list(f.factornf(R.pari_polynomial("a")))
     2880            # PARI's nffactor() ignores the unit, _factor_pari_helper()
     2881            # adds back the unit of the factorization.
     2882            return self._factor_pari_helper(G)
    28912883
    28922884        elif is_RealField(R):
    28932885            n = pari.set_real_precision(int(3.5*R.prec()) + 1)