Line data Source code
1 : /* Copyright (C) 2000-2004 The PARI group.
2 :
3 : This file is part of the PARI/GP package.
4 :
5 : PARI/GP is free software; you can redistribute it and/or modify it under the
6 : terms of the GNU General Public License as published by the Free Software
7 : Foundation; either version 2 of the License, or (at your option) any later
8 : version. It is distributed in the hope that it will be useful, but WITHOUT
9 : ANY WARRANTY WHATSOEVER.
10 :
11 : Check the License for details. You should have received a copy of it, along
12 : with the package; see the file 'COPYING'. If not, write to the Free Software
13 : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
14 :
15 : /***********************************************************************/
16 : /** **/
17 : /** ARITHMETIC OPERATIONS ON POLYNOMIALS **/
18 : /** (first part) **/
19 : /** **/
20 : /***********************************************************************/
21 : #include "pari.h"
22 : #include "paripriv.h"
23 : /*******************************************************************/
24 : /* */
25 : /* POLYNOMIAL EUCLIDEAN DIVISION */
26 : /* */
27 : /*******************************************************************/
28 : /* x t_POLMOD, y t_POL in the same variable as x[1], return x % y */
29 : static GEN
30 12887 : polmod_mod(GEN x, GEN y)
31 : {
32 12887 : GEN z, a, T = gel(x,1);
33 12887 : if (RgX_equal(T, y)) return gcopy(x);
34 14 : z = cgetg(3,t_POLMOD); T = RgX_gcd(T,y); a = gel(x,2);
35 14 : gel(z,1) = T;
36 14 : gel(z,2) = (typ(a)==t_POL && varn(a)==varn(T))? RgX_rem(a, T): gcopy(a);
37 14 : return z;
38 : }
39 : /* x,y two "scalars", return 0 with type info */
40 : static GEN
41 1575 : rem_scal_scal(GEN x, GEN y)
42 : {
43 1575 : pari_sp av = avma;
44 1575 : GEN z = gadd(gmul(gen_0,x), gmul(gen_0,y));
45 1575 : if (gequal0(y)) pari_err_INV("grem",y);
46 1575 : return gerepileupto(av, simplify(z));
47 : }
48 : /* x pol, y "scalar", return 0 with type info */
49 : static GEN
50 126 : rem_pol_scal(GEN x, GEN y)
51 : {
52 126 : pari_sp av = avma;
53 126 : if (gequal0(y)) pari_err_INV("grem",y);
54 126 : return gerepileupto(av, simplify(gmul(Rg_get_0(x),y)));
55 : }
56 : /* x "scalar", y pol, return x % y with type info */
57 : static GEN
58 1054908 : rem_scal_pol(GEN x, GEN y)
59 : {
60 1054908 : if (degpol(y))
61 : {
62 1053340 : if (!signe(y)) pari_err_INV("grem",y);
63 1053340 : return gmul(x, Rg_get_1(y));
64 : }
65 1568 : y = gel(y,2); return rem_scal_scal(x,y);
66 : }
67 : GEN
68 273 : poldivrem(GEN x, GEN y, GEN *pr)
69 : {
70 273 : const char *f = "euclidean division";
71 273 : long tx = typ(x), ty = typ(y), vx = gvar(x), vy = gvar(y);
72 : GEN z;
73 :
74 273 : if (!is_extscalar_t(tx) || !is_extscalar_t(ty)) pari_err_TYPE2(f,x,y);
75 273 : if (vx == vy && ((tx==t_POLMOD) ^ (ty==t_POLMOD))) pari_err_TYPE2(f,x,y);
76 259 : if (ty != t_POL || varncmp(vx, vy) < 0) /* y "scalar" */
77 : {
78 70 : if (!pr || pr == ONLY_DIVIDES) return gdiv(x,y);
79 70 : if (tx != t_POL || varncmp(vy, vx) < 0) /* x "scalar" */
80 0 : z = rem_scal_scal(x,y);
81 : else
82 70 : z = rem_pol_scal(x,y);
83 70 : if (pr == ONLY_REM) return z;
84 70 : *pr = z; return gdiv(x,y);
85 : }
86 189 : if (tx != t_POL || varncmp(vx, vy) > 0) /* x "scalar" */
87 : {
88 84 : if (!degpol(y)) /* constant t_POL, treat as scalar */
89 : {
90 7 : y = gel(y,2);
91 7 : if (!pr || pr == ONLY_DIVIDES) gdiv(x,y);
92 7 : z = rem_scal_scal(x,y);
93 7 : if (pr == ONLY_REM) return z;
94 7 : *pr = z; return gdiv(x,y);
95 : }
96 77 : if (!signe(y)) pari_err_INV("poldivrem",y);
97 77 : if (!pr || pr == ONLY_DIVIDES) return gequal0(x)? Rg_get_0(y): NULL;
98 77 : z = gmul(x, Rg_get_1(y));
99 77 : if (pr == ONLY_REM) return z;
100 77 : *pr = z; return Rg_get_0(y);
101 : }
102 105 : return RgX_divrem(x,y,pr);
103 : }
104 : GEN
105 637 : gdeuc(GEN x, GEN y)
106 : {
107 637 : const char *f = "euclidean division";
108 637 : long tx = typ(x), ty = typ(y), vx = gvar(x), vy = gvar(y);
109 637 : if (!is_extscalar_t(tx) || !is_extscalar_t(ty)) pari_err_TYPE2(f,x,y);
110 623 : if (vx == vy && ((tx==t_POLMOD) ^ (ty==t_POLMOD))) pari_err_TYPE2(f,x,y);
111 595 : if (ty != t_POL || varncmp(vx, vy) < 0) return gdiv(x,y); /* y "scalar" */
112 455 : if (tx != t_POL || varncmp(vx, vy) > 0)
113 : { /* x "scalar" */
114 140 : if (!signe(y)) pari_err_INV("gdeuc",y);
115 140 : if (!degpol(y)) return gdiv(x, gel(y,2)); /* constant */
116 140 : return Rg_get_0(y);
117 : }
118 315 : return RgX_div(x,y);
119 : }
120 : GEN
121 4199587 : grem(GEN x, GEN y)
122 : {
123 4199587 : const char *f = "euclidean division";
124 4199587 : long tx = typ(x), ty = typ(y), vx = gvar(x), vy = gvar(y);
125 :
126 4199584 : if (ty == t_POL)
127 : {
128 4199521 : if (varncmp(vx,vy) >= 0)
129 : {
130 : pari_sp av;
131 : GEN z;
132 4199521 : if (!signe(y)) pari_err_INV("grem",y);
133 4199524 : if (vx != vy) return rem_scal_pol(x,y);
134 3144616 : switch(tx)
135 : {
136 12887 : case t_POLMOD: return polmod_mod(x,y);
137 3119598 : case t_POL: return RgX_rem(x,y);
138 12082 : case t_RFRAC:
139 : {
140 12082 : GEN a = gel(x,1), b = gel(x,2), p, pol;
141 12082 : if (typ(a) == t_POL && RgX_is_ZX(y) && ZX_is_monic(y))
142 : {
143 12047 : long pa, t = RgX_type2(a,b, &p,&pol,&pa);
144 12047 : if (t == t_FRAC || t == t_INT) return QXQ_div(a, b, y);
145 : }
146 35 : av = avma; z = RgXQ_inv(RgX_rem(b, y), y);
147 28 : return gerepileupto(av, grem(gmul(a, z), y));
148 : }
149 49 : case t_SER:
150 49 : if (RgX_is_monomial(y))
151 : {
152 28 : if (lg(x)-2 + valser(x) < degpol(y)) pari_err_OP("%",x,y);
153 21 : av = avma;
154 21 : return gerepileupto(av, gmod(ser2rfrac_i(x), y));
155 : }
156 21 : default: pari_err_TYPE2("%",x,y);
157 : }
158 : }
159 0 : else switch(tx)
160 : {
161 0 : case t_POL:
162 0 : case t_RFRAC: return rem_pol_scal(x,y);
163 0 : default: pari_err_TYPE2("%",x,y);
164 : }
165 : }
166 63 : if (!is_extscalar_t(tx) || !is_extscalar_t(ty)) pari_err_TYPE2(f,x,y);
167 63 : if (vx == vy && ty==t_POLMOD) pari_err_TYPE2(f,x,y);
168 56 : if (tx != t_POL || varncmp(vx,vy) > 0)
169 : { /* x a "scalar" */
170 0 : if (ty != t_POL || varncmp(vx, vy) < 0) return rem_scal_scal(x,y);
171 0 : return rem_scal_pol(x,y);
172 : }
173 56 : if (ty != t_POL || varncmp(vx, vy) < 0) /* y a "scalar" */
174 56 : return rem_pol_scal(x,y);
175 0 : return RgX_rem(x,y);
176 : }
177 :
178 : /*******************************************************************/
179 : /* */
180 : /* CONVERSIONS RELATED TO p-ADICS */
181 : /* */
182 : /*******************************************************************/
183 : /* x t_PADIC, p a prime or NULL (unset). Consistency check */
184 : static void
185 336 : check_padic_p(GEN x, GEN p)
186 : {
187 336 : GEN q = gel(x,2);
188 336 : if (p && !equalii(p, q)) pari_err_MODULUS("Zp_to_Z", p,q);
189 315 : }
190 : /* shallow */
191 : static GEN
192 4326 : Zp_to_Z(GEN x, GEN p) {
193 4326 : switch(typ(x))
194 : {
195 4088 : case t_INT: break;
196 238 : case t_PADIC:
197 238 : check_padic_p(x, p);
198 217 : x = gtrunc(x); break;
199 0 : default: pari_err_TYPE("Zp_to_Z",x);
200 : }
201 4305 : return x;
202 : }
203 : /* shallow */
204 : static GEN
205 749 : ZpX_to_ZX(GEN x, GEN p)
206 4914 : { pari_APPLY_pol_normalized(Zp_to_Z(gel(x,i), p)); }
207 :
208 : static GEN
209 686 : get_padic_content(GEN f, GEN p)
210 : {
211 686 : GEN c = content(f);
212 686 : if (gequal0(c)) /* O(p^n) can occur */
213 : {
214 0 : if (typ(c) != t_PADIC) pari_err_TYPE("QpX_to_ZX",f);
215 0 : check_padic_p(c, p);
216 0 : c = powis(p, valp(c));
217 : }
218 686 : return c;
219 : }
220 : /* make f suitable for [root|factor]padic. Shallow */
221 : static GEN
222 623 : QpX_to_ZX(GEN f, GEN p)
223 : {
224 623 : GEN c = get_padic_content(f, p);
225 623 : f = RgX_Rg_div(f, c);
226 623 : return ZpX_to_ZX(f, p);
227 : }
228 :
229 : /* x in Z return x + O(pr), pr = p^r. Shallow */
230 : static GEN
231 4361 : Z_to_Zp(GEN x, GEN p, GEN pr, long r)
232 : {
233 : GEN y;
234 4361 : long v, sx = signe(x);
235 :
236 4361 : if (!sx) return zeropadic_shallow(p,r);
237 3829 : v = Z_pvalrem(x,p,&x);
238 3829 : if (v) {
239 840 : if (r <= v) return zeropadic_shallow(p,r);
240 735 : r -= v;
241 735 : pr = powiu(p,r);
242 : }
243 3724 : y = cgetg(5,t_PADIC);
244 3724 : y[1] = evalprecp(r)|evalvalp(v);
245 3724 : gel(y,2) = p;
246 3724 : gel(y,3) = pr;
247 3724 : gel(y,4) = modii(x,pr); return y;
248 : }
249 : /* shallow */
250 : static GEN
251 56 : ZV_to_ZpV(GEN z, GEN p, long r)
252 : {
253 56 : long i, l = lg(z);
254 56 : GEN Z = cgetg(l, typ(z)), q = powiu(p, r);
255 161 : for (i=1; i<l; i++) gel(Z,i) = Z_to_Zp(gel(z,i),p,q,r);
256 56 : return Z;
257 : }
258 : /* shallow */
259 : static GEN
260 1253 : ZX_to_ZpX(GEN z, GEN p, GEN q, long r)
261 : {
262 1253 : long i, l = lg(z);
263 1253 : GEN Z = cgetg(l, t_POL); Z[1] = z[1];
264 5327 : for (i=2; i<l; i++) gel(Z,i) = Z_to_Zp(gel(z,i),p,q,r);
265 1253 : return Z;
266 : }
267 : /* return (x + O(p^r)) normalized (multiply by a unit such that leading coeff
268 : * is a power of p), x in Z[X] (or Z_p[X]). Shallow */
269 : static GEN
270 1169 : ZX_to_ZpX_normalized(GEN x, GEN p, GEN pr, long r)
271 : {
272 1169 : long i, lx = lg(x);
273 1169 : GEN z, lead = leading_coeff(x);
274 :
275 1169 : if (gequal1(lead)) return ZX_to_ZpX(x, p, pr, r);
276 56 : (void)Z_pvalrem(lead, p, &lead); lead = Fp_inv(lead, pr);
277 56 : z = cgetg(lx,t_POL);
278 238 : for (i=2; i < lx; i++) gel(z,i) = Z_to_Zp(mulii(lead,gel(x,i)),p,pr,r);
279 56 : z[1] = x[1]; return z;
280 : }
281 : static GEN
282 49 : ZXV_to_ZpXQV(GEN z, GEN T, GEN p, long r)
283 : {
284 49 : long i, l = lg(z);
285 49 : GEN Z = cgetg(l, typ(z)), q = powiu(p, r);
286 49 : T = ZX_copy(T);
287 126 : for (i=1; i<lg(z); i++) gel(Z,i) = mkpolmod(ZX_to_ZpX(gel(z,i),p,q,r),T);
288 49 : return Z;
289 : }
290 : /* shallow */
291 : static GEN
292 63 : QpXQX_to_ZXY(GEN f, GEN p)
293 : {
294 63 : GEN c = get_padic_content(f, p);
295 63 : long i, l = lg(f);
296 63 : f = RgX_Rg_div(f,c);
297 287 : for (i=2; i<l; i++)
298 : {
299 231 : GEN t = gel(f,i);
300 231 : switch(typ(t))
301 : {
302 91 : case t_POLMOD:
303 91 : t = gel(t,2);
304 91 : t = (typ(t) == t_POL)? ZpX_to_ZX(t, p): Zp_to_Z(t, p);
305 91 : break;
306 0 : case t_POL: t = ZpX_to_ZX(t, p); break;
307 140 : default: t = Zp_to_Z(t, p); break;
308 : }
309 224 : gel(f,i) = t;
310 : }
311 56 : return f;
312 : }
313 :
314 : /*******************************************************************/
315 : /* */
316 : /* p-ADIC ROOTS */
317 : /* */
318 : /*******************************************************************/
319 :
320 : /* f primitive ZX, squarefree, leading term prime to p; 0 <= a < p such that
321 : * f(a) = 0 mod p. Return p-adic roots of f equal to a mod p, in
322 : * precision >= prec */
323 : GEN
324 2863 : ZX_Zp_root(GEN f, GEN a, GEN p, long prec)
325 : {
326 : GEN z, R;
327 : long i, j, k;
328 :
329 2863 : if (signe(FpX_eval(FpX_deriv(f, p), a, p)))
330 : { /* simple zero mod p, go all the way to p^prec */
331 2632 : if (prec > 1) a = ZpX_liftroot(f, a, p, prec);
332 2632 : return mkcol(a);
333 : }
334 :
335 231 : f = ZX_unscale_div(ZX_translate(f,a), p); /* f(pX + a) / p */
336 231 : (void)ZX_pvalrem(f,p,&f);
337 231 : z = cgetg(degpol(f)+1,t_COL);
338 :
339 231 : R = FpX_roots(f, p);
340 574 : for (j=i=1; i<lg(R); i++)
341 : {
342 343 : GEN u = ZX_Zp_root(f, gel(R,i), p, prec-1);
343 756 : for (k=1; k<lg(u); k++) gel(z,j++) = addii(a, mulii(p, gel(u,k)));
344 : }
345 231 : setlg(z,j); return z;
346 : }
347 :
348 : /* a t_PADIC, return vector of p-adic roots of f equal to a (mod p) */
349 : GEN
350 56 : Zp_appr(GEN f, GEN a)
351 : {
352 56 : pari_sp av = avma;
353 56 : GEN z, p = gel(a,2);
354 56 : long v = valp(a), prec = v;
355 :
356 56 : if (signe(gel(a,4))) prec += precp(a);
357 56 : f = QpX_to_ZX(f, p);
358 42 : if (degpol(f) <= 0) pari_err_CONSTPOL("Zp_appr");
359 42 : if (v < 0) pari_err_DOMAIN("padicappr", "v(a)", "<", gen_0, stoi(v));
360 35 : f = ZX_radical(f);
361 35 : a = padic_to_Fp(a, p);
362 35 : if (signe(FpX_eval(f,a,p))) { set_avma(av); return cgetg(1,t_COL); }
363 28 : z = ZX_Zp_root(f, a, p, prec);
364 28 : return gerepilecopy(av, ZV_to_ZpV(z, p, prec));
365 : }
366 : static long
367 126 : pval(GEN x, GEN p) { return typ(x) == t_INT? Z_pval(x,p): ZX_pval(x,p); }
368 : /* f a ZXX, f(0) != 0 */
369 : static GEN
370 539 : pnormalize(GEN f, GEN p, GEN T, long prec, long n,
371 : GEN *plead, long *pprec, int *prev)
372 : {
373 539 : *plead = leading_coeff(f);
374 539 : *pprec = prec;
375 539 : *prev = 0;
376 539 : if (!isint1(*plead))
377 : {
378 63 : long v = pval(*plead,p), v1 = pval(constant_coeff(f),p);
379 63 : if (v1 < v)
380 : {
381 42 : *prev = 1;
382 42 : f = RgX_recip_i(f); /* f(0) != 0 so degree is the same */
383 : /* beware loss of precision from lc(factor), whose valuation is <= v */
384 42 : *pprec += v; v = v1;
385 : }
386 63 : *pprec += v * n;
387 : }
388 539 : if (!T) return ZX_Q_normalize(f, plead);
389 14 : *plead = gen_1;
390 14 : return FpXQX_normalize(f, T, powiu(p,*pprec));
391 : }
392 :
393 : /**************************************************************************/
394 :
395 : static void
396 238 : scalar_getprec(GEN x, long *pprec, GEN *pp)
397 : {
398 238 : if (typ(x)==t_PADIC)
399 : {
400 98 : long e = valp(x); if (signe(gel(x,4))) e += precp(x);
401 98 : if (e < *pprec) *pprec = e;
402 98 : check_padic_p(x, *pp);
403 98 : *pp = gel(x,2);
404 : }
405 238 : }
406 : static void
407 98 : getprec(GEN x, long *pprec, GEN *pp)
408 : {
409 : long i;
410 98 : if (typ(x) != t_POL) scalar_getprec(x, pprec, pp);
411 : else
412 266 : for (i = lg(x)-1; i>1; i--) scalar_getprec(gel(x,i), pprec, pp);
413 98 : }
414 :
415 : /* assume f(a) = 0 (mod T,p) */
416 : static GEN
417 105 : ZXY_ZpQ_root(GEN f, GEN a, GEN T, GEN p, long prec)
418 : {
419 : GEN z, R;
420 : long i, j, k, lR;
421 105 : if (signe(FqX_eval(FqX_deriv(f,T,p), a, T,p)))
422 : { /* simple zero mod (T,p), go all the way to p^prec */
423 77 : if (prec > 1) a = ZpXQX_liftroot(f, a, T, p, prec);
424 77 : return mkcol(a);
425 : }
426 28 : f = RgX_unscale(RgXQX_translate(f, a, T), p);
427 28 : f = RgX_Rg_div(f, powiu(p, gvaluation(f,p)));
428 28 : z = cgetg(degpol(f)+1,t_COL);
429 28 : R = FpXQX_roots(FqX_red(f,T,p), T, p); lR = lg(R);
430 70 : for(j=i=1; i<lR; i++)
431 : {
432 42 : GEN u = ZXY_ZpQ_root(f, gel(R,i), T, p, prec-1);
433 84 : for (k=1; k<lg(u); k++) gel(z,j++) = gadd(a, gmul(p, gel(u,k)));
434 : }
435 28 : setlg(z,j); return z;
436 : }
437 :
438 : /* a belongs to finite extension of Q_p, return all roots of f equal to a
439 : * mod p. Don't assume f(a) = 0 (mod p) */
440 : GEN
441 105 : padicappr(GEN f, GEN a)
442 : {
443 : GEN p, z, T, Tp;
444 : long prec;
445 105 : pari_sp av = avma;
446 :
447 105 : if (typ(f)!=t_POL) pari_err_TYPE("padicappr",f);
448 105 : switch(typ(a)) {
449 56 : case t_PADIC: return Zp_appr(f,a);
450 49 : case t_POLMOD: break;
451 0 : default: pari_err_TYPE("padicappr",a);
452 : }
453 49 : if (gequal0(f)) pari_err_ROOTS0("padicappr");
454 49 : T = gel(a,1);
455 49 : a = gel(a,2);
456 49 : p = NULL; prec = LONG_MAX;
457 49 : getprec(a, &prec, &p);
458 49 : getprec(T, &prec, &p); if (!p) pari_err_TYPE("padicappr",T);
459 49 : f = QpXQX_to_ZXY(f, p);
460 42 : if (typ(a) != t_POL) a = scalarpol_shallow(a, varn(T));
461 42 : a = ZpX_to_ZX(a,p);
462 42 : T = QpX_to_ZX(T,p);
463 : /* ensure that f /= (f,f') is separable */
464 42 : (void)nfgcd_all(f, RgX_deriv(f), T, NULL, &f);
465 :
466 42 : Tp = FpX_red(T, p); a = FqX_red(a, Tp, p);
467 42 : if (!gequal0(FqX_eval(FqX_red(f,Tp,p), a, Tp,p)))
468 7 : { set_avma(av); return cgetg(1,t_COL); } /* f(a) != 0 (mod p,T) */
469 35 : z = ZXY_ZpQ_root(f, a, T, p, prec);
470 35 : return gerepilecopy(av, ZXV_to_ZpXQV(z, T, p, prec));
471 : }
472 :
473 : /* vector of p-adic roots of the ZX f, leading term prime to p. Shallow */
474 : static GEN
475 35 : ZX_Zp_roots(GEN f, GEN p, long prec)
476 : {
477 : long l, i;
478 : GEN r;
479 :
480 35 : f = ZX_radical(f);
481 35 : r = FpX_roots(f, p);
482 35 : l = lg(r); if (l == 1) return r;
483 91 : for (i = 1; i < l; i++) gel(r,i) = ZX_Zp_root(f, gel(r,i), p, prec);
484 28 : return ZV_to_ZpV(shallowconcat1(r), p, prec);
485 : }
486 : /* vector of p-adic roots of the ZXX f, leading term prime to p. Shallow */
487 : static GEN
488 14 : ZXY_ZpQ_roots(GEN f, GEN T, GEN p, long prec)
489 : {
490 : long l, i;
491 : GEN r;
492 :
493 14 : (void)nfgcd_all(f, RgX_deriv(f), T, NULL, &f);
494 14 : r = FqX_roots(f, FpX_red(T,p), p);
495 14 : l = lg(r); if (l == 1) return r;
496 42 : for (i = 1; i < l; i++) gel(r,i) = ZXY_ZpQ_root(f, gel(r,i), T, p, prec);
497 14 : return ZXV_to_ZpXQV(shallowconcat1(r), T, p, prec);
498 : }
499 :
500 : /* return p-adic roots of f, precision prec */
501 : GEN
502 56 : polrootspadic(GEN f, GEN Tp, long prec)
503 : {
504 56 : pari_sp av = avma;
505 : GEN lead, y, T, p;
506 : long PREC, i, k, v;
507 : int reverse;
508 :
509 56 : if (!ff_parse_Tp(Tp, &T,&p,0)) pari_err_TYPE("polrootspadic",Tp);
510 56 : if (typ(f)!=t_POL) pari_err_TYPE("polrootspadic",f);
511 56 : if (gequal0(f)) pari_err_ROOTS0("polrootspadic");
512 56 : if (prec <= 0)
513 7 : pari_err_DOMAIN("polrootspadic", "precision", "<=",gen_0,stoi(prec));
514 49 : f = T? QpXQX_to_ZXY(f, p): QpX_to_ZX(f, p);
515 49 : v = RgX_valrem(f, &f);
516 49 : f = pnormalize(f, p, T, prec, 1, &lead, &PREC, &reverse);
517 49 : y = T? ZXY_ZpQ_roots(f,T,p,PREC): ZX_Zp_roots(f,p,PREC);
518 49 : k = lg(y);
519 49 : if (lead != gen_1) y = RgC_Rg_div(y, lead);
520 49 : if (reverse)
521 7 : for (i=1; i<k; i++) gel(y,i) = ginv(gel(y,i));
522 49 : if (v) y = shallowconcat(zeropadic_shallow(p, prec), y);
523 49 : return gerepilecopy(av, y);
524 : }
525 :
526 : /*******************************************************************/
527 : /* */
528 : /* FACTORIZATION in Zp[X], using ROUND4 */
529 : /* */
530 : /*******************************************************************/
531 :
532 : int
533 3007 : cmp_padic(GEN x, GEN y)
534 : {
535 : long vx, vy;
536 3007 : if (x == gen_0) return -1;
537 3007 : if (y == gen_0) return 1;
538 3007 : vx = valp(x);
539 3007 : vy = valp(y);
540 3007 : if (vx < vy) return 1;
541 2972 : if (vx > vy) return -1;
542 2734 : return cmpii(gel(x,4), gel(y,4));
543 : }
544 :
545 : /* replace p^e by p*...*p [ factors are not known to be equal, only close at
546 : * input accuracy ] */
547 : static GEN
548 49 : famat_flatten(GEN fa)
549 : {
550 49 : GEN y, P = gel(fa,1), E = gel(fa,2);
551 49 : long i, l = lg(E);
552 49 : y = cgetg(l, t_VEC);
553 147 : for (i = 1; i < l; i++)
554 : {
555 98 : GEN p = gel(P,i);
556 98 : long e = itou(gel(E,i));
557 98 : gel(y,i) = const_col(e, p);
558 : }
559 49 : y = shallowconcat1(y); return mkmat2(y, const_col(lg(y)-1, gen_1));
560 : }
561 :
562 : GEN
563 525 : factorpadic(GEN f, GEN p, long r)
564 : {
565 525 : pari_sp av = avma;
566 : GEN y, ppow;
567 : long v, n;
568 525 : int reverse = 0, exact;
569 :
570 525 : if (typ(f)!=t_POL) pari_err_TYPE("factorpadic",f);
571 525 : if (typ(p)!=t_INT) pari_err_TYPE("factorpadic",p);
572 525 : if (r <= 0) pari_err_DOMAIN("factorpadic", "precision", "<=",gen_0,stoi(r));
573 518 : if (!signe(f)) return prime_fact(f);
574 518 : if (!degpol(f)) return trivial_fact();
575 :
576 518 : exact = RgX_is_QX(f); /* before RgX_valrem which may lose type information */
577 518 : v = RgX_valrem_inexact(f, &f);
578 518 : ppow = powiu(p,r);
579 518 : n = degpol(f);
580 518 : if (!n)
581 28 : y = trivial_fact();
582 : else
583 : {
584 : GEN P, lead;
585 : long i, l, pr;
586 :
587 490 : f = QpX_to_ZX(f, p);
588 490 : f = pnormalize(f, p, NULL, r, n-1, &lead, &pr, &reverse);
589 490 : y = ZpX_monic_factor(f, p, pr);
590 490 : P = gel(y,1); l = lg(P);
591 490 : if (!isint1(lead))
592 266 : for (i=1; i<l; i++) gel(P,i) = Q_primpart(RgX_unscale(gel(P,i), lead));
593 1659 : for (i=1; i<l; i++)
594 : {
595 1169 : GEN t = gel(P,i);
596 1169 : if (reverse) t = RgX_recip_shallow(t);
597 1169 : gel(P,i) = ZX_to_ZpX_normalized(t,p,ppow,r);
598 : }
599 : }
600 518 : if (v)
601 : { /* v > 0 */
602 63 : GEN X = ZX_to_ZpX(pol_x(varn(f)), p, ppow, r);
603 63 : y = famat_mulpow_shallow(y, X, utoipos(v));
604 : }
605 518 : if (!exact) y = famat_flatten(y);
606 518 : return gerepilecopy(av, sort_factor_pol(y, cmp_padic));
607 : }
|
