Ticket #10910: factor_nfinit_free.patch

File factor_nfinit_free.patch, 5.9 KB (added by lftabera, 9 years ago)

the good one

  • sage/rings/polynomial/polynomial_element.pyx

    # HG changeset patch
    # User Luis Felipe Tabera Alonso <lftabera@yahoo.es>
    # Parent 2a2abbcad325ccca9399981ceddf5897eb467e64
    #10910 Avoid nfinit while factoring polynomials
    
    diff -r 2a2abbcad325 sage/rings/polynomial/polynomial_element.pyx
    a b  
    27092709            sage: pol.factor()
    27102710            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)
    27112711
     2712        Check that we do not compute nfinit::
     2713
     2714            sage: x = var('x')
     2715            sage: N = NumberField(x^2-2,'x')
     2716            sage: K = N['x']
     2717            sage: f1 = K.random_element(3)
     2718            sage: f2 = K.random_element(3)
     2719            sage: f = f1*f2
     2720            sage: hasattr(N, '_NumberField_generic__pari_nf')
     2721            False
     2722            sage: F1 = factor(f)
     2723            sage: hasattr(N, '_NumberField_generic__pari_nf')
     2724            False
     2725            sage: _ = N.pari_nf()
     2726            sage: F2 = factor(f)
     2727            sage: hasattr(N, '_NumberField_generic__pari_nf')
     2728            True
     2729            sage: F1 -f
     2730            0
     2731            sage: F2 -f
     2732            0
     2733
     2734        Next example would be very slow if we had to compute a nfinit             structure. See #10910::
     2735
     2736            sage: x = var('x')
     2737            sage: N = NumberField(x^8+1212/413*x^7+1222*x+1/17,'a')
     2738            sage: K.<t> = N[]
     2739            sage: f = K.random_element(4)
     2740            sage: g = K.random_element(4)
     2741            sage: F = factor(f*g)
     2742            sage: F - f*g
     2743            0
     2744            sage: hasattr(N, '_NumberField_generic__pari_nf')
     2745            False
     2746
     2747        Check that the factorization is correct. See #10910::
     2748
     2749            sage: K=QQ[sqrt(2**40+1), 3^(1/3)]
     2750            sage: x = polygen(K)
     2751            sage: f = x^3 - 3
     2752            sage: len(factor(f))
     2753            2
    27122754        """
    27132755        # PERFORMANCE NOTE:
    27142756        #     In many tests with SMALL degree PARI is substantially
     
    27772819            g = M['x']([to_M(x) for x in self.list()])
    27782820            F = g.factor()
    27792821            S = self.parent()
    2780             v = [(S([from_M(x) for x in f.list()]), e) for f, e in g.factor()]
     2822            v = [(S([from_M(x) for x in f.list()]), e) for f, e in F]
    27812823            return Factorization(v, from_M(F.unit()))
    27822824
    27832825        elif is_FiniteField(R):
     
    28112853            if R.defining_polynomial().denominator() == 1:
    28122854                v = [ x._pari_("a") for x in self.list() ]
    28132855                f = pari(v).Polrev()
    2814                 Rpari = R.pari_nf()
    2815                 if (Rpari.variable() != "a"):
    2816                     Rpari = copy.copy(Rpari)
    2817                     Rpari[0] = Rpari[0]("a")
    2818                     Rpari[6] = [ x("a") for x in Rpari[6] ]
    2819                 G = list(Rpari.nffactor(f))
    2820                 # PARI's nffactor() ignores the unit, _factor_pari_helper()
    2821                 # adds back the unit of the factorization.
     2856                if hasattr(R, '_NumberField_generic__pari_nf'):
     2857                    # If we have a nf structure, this should be faster
     2858                    Rpari = R.pari_nf()
     2859                    if (Rpari.variable() != "a"):
     2860                        Rpari = copy.copy(Rpari)
     2861                        Rpari[0] = Rpari[0]("a")
     2862                        Rpari[6] = [ x("a") for x in Rpari[6] ]
     2863                    G = list(Rpari.nffactor(f))
     2864                    # PARI's nffactor() ignores the unit, _factor_pari_helper()
     2865                    # adds back the unit of the factorization.
     2866
     2867                else:
     2868                    # We do not compute a nf structure, since it may be
     2869                    # expensive
     2870                    nf = R.pari_polynomial('a')
     2871                    G = list(nf.nffactor(f))
     2872
    28222873                return self._factor_pari_helper(G)
    28232874
    28242875            else:
    28252876
    28262877                Rdenom = R.defining_polynomial().denominator()
    2827 
    2828                 new_Rpoly = (R.defining_polynomial() * Rdenom).change_variable_name("a")
    2829 
    2830                 Rpari, Rdiff = new_Rpoly._pari_().nfinit(3)
    2831 
    2832                 AZ = QQ['z']
    2833                 Raux = NumberField(AZ(Rpari[0]),'alpha')
    2834 
    2835                 S, gSRaux, fRauxS = Raux.change_generator(Raux(Rdiff))
    2836 
    2837                 phi_RS = R.Hom(S)([S.gen(0)])
    2838                 phi_SR = S.Hom(R)([R.gen(0)])
    2839 
     2878                d = R.degree()
     2879                Rdenom_d = Rdenom**d
     2880                x = R.polynomial_ring().gen()
     2881                new_Rpoly = R.defining_polynomial()(x/Rdenom)*Rdenom_d
     2882
     2883                S = NumberField(new_Rpoly, 'alpha')
     2884                phi_RS = R.Hom(S)([S.gen()/Rdenom])
     2885                phi_SR = S.Hom(R)([R.gen()*Rdenom])
    28402886                unit = self.leading_coefficient()
    2841                 temp_f = self * 1/unit
    2842 
    2843                 v = [ gSRaux(phi_RS(x))._pari_("a") for x in temp_f.list() ]
     2887                temp_f = self * ~unit
     2888                v = [ (phi_RS(x))._pari_("a") for x in temp_f.list() ]
    28442889                f = pari(v).Polrev()
    2845 
    2846                 pari_factors = Rpari.nffactor(f)
    2847 
    2848                 factors = [ ( self.parent([ phi_SR(fRauxS(Raux(pari_factors[0][i][j])))
    2849                                             for j in range(len(pari_factors[0][i])) ]) ,
    2850                              int(pari_factors[1][i]) )
     2890                nf = S.pari_polynomial('a')
     2891                pari_factors = nf.nffactor(f)
     2892                factors = [ ( self.parent([ phi_SR(S(pari_factors[0][i][j]))
     2893                            for j in range(len(pari_factors[0][i])) ]) ,
     2894                            int(pari_factors[1][i]) )
    28512895                            for i in range(pari_factors.nrows()) ]
    28522896
    28532897                return Factorization(factors, unit)
    2854        
    28552898
    28562899        elif is_RealField(R):
    28572900            n = pari.set_real_precision(int(3.5*R.prec()) + 1)