1"""Square-free decomposition algorithms and related tools."""
2
3from .polyerrors import DomainError
4
5
6class _SQF:
7    """Mixin class to compute square-free decomposition of polynomials."""
8
9    def sqf_list(self, f):
10        """
11        Return square-free decomposition of a polynomial in ``K[X]``.
12
13        Examples
14        ========
15
16        >>> R, x, y = ring('x y', ZZ)
17
18        >>> R.sqf_list(x**5 + 2*x**4*y + x**3*y**2)
19        (1, [(x + y, 2), (x, 3)])
20
21        >>> R, x, y = ring('x y', FF(5))
22
23        >>> f = x**5*y**5 + 1
24
25        Note that e.g.:
26
27        >>> f.diff(x)
28        0 mod 5
29
30        This phenomenon doesn't happen in characteristic zero. However we can
31        still compute square-free decomposition of ``f``:
32
33        >>> R.sqf_list(f)
34        (1 mod 5, [(x*y + 1 mod 5, 5)])
35
36        """
37        domain = self.domain
38
39        if domain.is_Field:
40            coeff, f = f.LC, f.monic()
41        else:
42            coeff, f = f.primitive()
43
44        if domain.is_FiniteField:
45            return coeff, self._gf_sqf_list(f)
46        else:
47            return coeff, self._rr_yun0_sqf_list(f)
48
49    def _gf_sqf_list(self, f):
50        """
51        Compute square-free decomposition of the monic ``f`` in ``GF(q)[X]``.
52
53        Notes
54        =====
55
56        Uses a modified version of Musser's algorithm for square-free
57        decomposition of univariate polynomials over finite fields.
58
59        References
60        ==========
61
62        * :cite:`Geddes1992algorithms`, algorithm 8.3
63        * :cite:`Musser1971factor`, p.54, algorithm V (modified)
64
65        """
66        domain = self.domain
67
68        n, factors, p = 1, [], int(domain.characteristic)
69        m = int(domain.order // p)
70
71        while not f.is_ground:
72            df = [f.diff(x) for x in self.gens]
73
74            if any(_ for _ in df):
75                g = f
76                for q in df:
77                    g = self.gcd(g, q)
78                h, f, i = f // g, g, 1
79
80                while h != 1:
81                    g = self.gcd(f, h)
82                    h //= g
83
84                    if not h.is_ground:
85                        factors.append((h, i*n))
86
87                    f //= g
88                    h = g
89                    i += 1
90
91            n *= p
92
93            g = self.zero
94            for monom, coeff in f.items():
95                g[(_ // p for _ in monom)] = coeff**m
96            f = g
97
98        return factors
99
100    def _rr_yun0_sqf_list(self, f):
101        """Compute square-free decomposition of ``f`` in zero-characteristic ring ``K[X]``.
102
103        References
104        ==========
105
106        * :cite:`Geddes1992algorithms`, algorithm 8.2
107        * :cite:`LeeM2013factor`, algorithm 3.1
108
109        """
110        if f.is_ground:
111            return []
112
113        result, count = [], 1
114        qs = [f.diff(x) for x in self.gens]
115
116        g = f
117        for q in qs:
118            g = self.gcd(g, q)
119
120        while f != 1:
121            qs = [q // g for q in qs]
122            f //= g
123            qs = [q - f.diff(x) for x, q in zip(self.gens, qs)]
124
125            g = f
126            for q in qs:
127                g = self.gcd(g, q)
128            if g != 1:
129                result.append((g, count))
130
131            count += 1
132
133        return result
134
135    def is_squarefree(self, f):
136        """
137        Return ``True`` if ``f`` is a square-free polynomial in ``K[X]``.
138
139        Examples
140        ========
141
142        >>> _, x, y = ring('x y', ZZ)
143
144        >>> ((x + y)**2).is_squarefree
145        False
146        >>> (x**2 + y**2).is_squarefree
147        True
148
149        """
150        if f.is_ground:
151            return True
152        else:
153            g = f
154            for x in self.gens:
155                g = self.gcd(g, f.diff(x))
156                if g.is_ground:
157                    return True
158            return False
159
160    def sqf_part(self, f):
161        """
162        Returns square-free part of a polynomial in ``K[X]``.
163
164        Examples
165        ========
166
167        >>> R, x, y = ring('x y', ZZ)
168
169        >>> R.sqf_part(x**3 + 2*x**2*y + x*y**2)
170        x**2 + x*y
171
172        """
173        domain = self.domain
174
175        if domain.is_FiniteField:
176            g = self.one
177            for f, _ in self.sqf_list(f)[1]:
178                g *= f
179
180            return g
181
182        if not f:
183            return f
184
185        gcd = f
186        for x in self.gens:
187            gcd = self.gcd(gcd, f.diff(x))
188        sqf = f // gcd
189
190        if domain.is_Field:
191            return sqf.monic()
192        else:
193            return sqf.primitive()[1]
194
195    def sqf_norm(self, f):
196        """
197        Square-free norm of ``f`` in ``K[X]``, useful over algebraic domains.
198
199        Returns ``s``, ``f``, ``r``, such that ``g(x) = f(x-sa)`` and ``r(x) = Norm(g(x))``
200        is a square-free polynomial over K, where ``a`` is the algebraic extension of ``K``.
201
202        Examples
203        ========
204
205        >>> _, x, y = ring('x y', QQ.algebraic_field(I))
206
207        >>> (x*y + y**2).sqf_norm()
208        (1, x*y - I*x + y**2 - 3*I*y - 2,
209         x**2*y**2 + x**2 + 2*x*y**3 + 2*x*y + y**4 + 5*y**2 + 4)
210
211        """
212        domain = self.domain
213
214        if not domain.is_AlgebraicField:
215            raise DomainError(f'ground domain must be algebraic, got {domain}')
216
217        new_ring = self.to_ground().inject(*domain.symbols, front=True)
218        g = domain.mod.set_ring(new_ring)
219        s = 0
220
221        while True:
222            h = f.inject(front=True)
223            r = g.resultant(h)
224
225            if r.is_squarefree:
226                return s, f, r
227            else:
228                f = f.compose({x: x - domain.unit for x in self.gens})
229                s += 1
230