Line data Source code
1 : /* Copyright (C) 2000 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 : /** GENERIC OPERATIONS **/
18 : /** (second part) **/
19 : /** **/
20 : /********************************************************************/
21 : #include "pari.h"
22 : #include "paripriv.h"
23 :
24 : /*********************************************************************/
25 : /** **/
26 : /** MAP FUNCTIONS WITH GIVEN PROTOTYPES **/
27 : /** **/
28 : /*********************************************************************/
29 : GEN
30 462 : map_proto_G(GEN (*f)(GEN), GEN x)
31 : {
32 462 : if (is_matvec_t(typ(x))) pari_APPLY_same(map_proto_G(f, gel(x,i)));
33 462 : return f(x);
34 : }
35 :
36 : GEN
37 13919578 : map_proto_lG(long (*f)(GEN), GEN x)
38 : {
39 13919662 : if (is_matvec_t(typ(x))) pari_APPLY_same(map_proto_lG(f, gel(x,i)));
40 13919550 : return stoi(f(x));
41 : }
42 :
43 : GEN
44 126 : map_proto_lGL(long (*f)(GEN,long), GEN x, long y)
45 : {
46 154 : if (is_matvec_t(typ(x))) pari_APPLY_same(map_proto_lGL(f,gel(x,i),y));
47 119 : return stoi(f(x,y));
48 : }
49 :
50 : static GEN
51 2258886 : _domul(void *data, GEN x, GEN y)
52 : {
53 2258886 : GEN (*mul)(GEN,GEN)=(GEN (*)(GEN,GEN)) data;
54 2258886 : return mul(x,y);
55 : }
56 :
57 : GEN
58 2462814 : gassoc_proto(GEN (*f)(GEN,GEN), GEN x, GEN y)
59 : {
60 2462814 : if (!y)
61 : {
62 2462814 : pari_sp av = avma;
63 2462814 : switch(typ(x))
64 : {
65 21 : case t_LIST:
66 21 : x = list_data(x); if (!x) return gen_1;
67 : case t_VEC:
68 2462800 : case t_COL: break;
69 7 : default: pari_err_TYPE("association",x);
70 : }
71 2462800 : return gc_upto(av, gen_product(x, (void *)f, _domul));
72 :
73 : }
74 0 : return f(x,y);
75 : }
76 :
77 : /*******************************************************************/
78 : /* */
79 : /* SIZES */
80 : /* */
81 : /*******************************************************************/
82 :
83 : long
84 5150313 : glength(GEN x)
85 : {
86 5150313 : long tx = typ(x);
87 5150313 : switch(tx)
88 : {
89 126 : case t_INT: return lgefint(x)-2;
90 609 : case t_LIST: {
91 609 : GEN L = list_data(x);
92 609 : return L? lg(L)-1: 0;
93 : }
94 14 : case t_REAL: return signe(x)? lg(x)-2: 0;
95 14 : case t_STR: return strlen( GSTR(x) );
96 91 : case t_VECSMALL: return lg(x)-1;
97 : }
98 5149459 : return lg(x) - lontyp[tx];
99 : }
100 :
101 : long
102 3878 : gtranslength(GEN x)
103 : {
104 3878 : switch(typ(x))
105 : {
106 0 : case t_VEC: case t_COL:
107 0 : return lg(x)-1;
108 3878 : case t_MAT:
109 3878 : return lg(x)==1 ? 0: nbrows(x);
110 0 : default:
111 0 : pari_err_TYPE("trans",x);
112 : return 0; /* LCOV_EXCL_LINE */
113 : }
114 : }
115 :
116 : GEN
117 1869 : matsize(GEN x)
118 : {
119 1869 : long L = lg(x) - 1;
120 1869 : switch(typ(x))
121 : {
122 7 : case t_VEC: return mkvec2s(1, L);
123 7 : case t_COL: return mkvec2s(L, 1);
124 1848 : case t_MAT: return mkvec2s(L? nbrows(x): 0, L);
125 : }
126 7 : pari_err_TYPE("matsize",x);
127 : return NULL; /* LCOV_EXCL_LINE */
128 : }
129 :
130 : /*******************************************************************/
131 : /* */
132 : /* CONVERSION GEN --> long */
133 : /* */
134 : /*******************************************************************/
135 :
136 : long
137 77 : gtolong(GEN x)
138 : {
139 77 : switch(typ(x))
140 : {
141 42 : case t_INT:
142 42 : return itos(x);
143 7 : case t_REAL:
144 7 : return (long)(rtodbl(x) + 0.5);
145 7 : case t_FRAC:
146 7 : { pari_sp av = avma; return gc_long(av, itos(ground(x))); }
147 7 : case t_COMPLEX:
148 7 : if (gequal0(gel(x,2))) return gtolong(gel(x,1)); break;
149 7 : case t_QUAD:
150 7 : if (gequal0(gel(x,3))) return gtolong(gel(x,2)); break;
151 : }
152 7 : pari_err_TYPE("gtolong",x);
153 : return 0; /* LCOV_EXCL_LINE */
154 : }
155 :
156 : /*******************************************************************/
157 : /* */
158 : /* COMPARISONS */
159 : /* */
160 : /*******************************************************************/
161 : static void
162 189 : chk_true_err()
163 : {
164 189 : GEN E = pari_err_last();
165 189 : switch(err_get_num(E))
166 : {
167 0 : case e_STACK: case e_MEM: case e_ALARM:
168 0 : pari_err(0, E); /* rethrow */
169 : }
170 189 : }
171 : /* x - y == 0 or undefined */
172 : static int
173 3195604 : gequal_try(GEN x, GEN y)
174 : {
175 : int i;
176 3195604 : pari_CATCH(CATCH_ALL) { chk_true_err(); return 0; }
177 3195604 : pari_TRY { i = gequal0(gadd(x, gneg_i(y))); } pari_ENDCATCH;
178 3195422 : return i;
179 : }
180 : /* x + y == 0 or undefined */
181 : static int
182 28 : gmequal_try(GEN x, GEN y)
183 : {
184 : int i;
185 28 : pari_CATCH(CATCH_ALL) { chk_true_err(); return 0; }
186 28 : pari_TRY { i = gequal0(gadd(x, y)); } pari_ENDCATCH;
187 21 : return i;
188 : }
189 :
190 : int
191 582925496 : isexactzero(GEN g)
192 : {
193 : long i, lx;
194 582925496 : switch (typ(g))
195 : {
196 477021609 : case t_INT:
197 477021609 : return !signe(g);
198 981974 : case t_INTMOD:
199 981974 : return !signe(gel(g,2));
200 15455448 : case t_COMPLEX:
201 15455448 : return isexactzero(gel(g,1)) && isexactzero(gel(g,2));
202 8111499 : case t_FFELT:
203 8111499 : return FF_equal0(g);
204 490 : case t_QUAD:
205 490 : return isexactzero(gel(g,2)) && isexactzero(gel(g,3));
206 305329 : case t_POLMOD:
207 305329 : return isexactzero(gel(g,2));
208 35578557 : case t_POL:
209 35578557 : lx = lg(g); /* cater for Mod(0,2)*x^0 */
210 35578557 : return lx == 2 || (lx == 3 && isexactzero(gel(g,2)));
211 461276 : case t_RFRAC:
212 461276 : return isexactzero(gel(g,1)); /* may occur: Mod(0,2)/x */
213 63504 : case t_VEC: case t_COL: case t_MAT:
214 63805 : for (i=lg(g)-1; i; i--)
215 63700 : if (!isexactzero(gel(g,i))) return 0;
216 105 : return 1;
217 : }
218 44945810 : return 0;
219 : }
220 : GEN
221 57171911 : gisexactzero(GEN g)
222 : {
223 : long i, lx;
224 : GEN a, b;
225 57171911 : switch (typ(g))
226 : {
227 24275684 : case t_INT:
228 24275684 : return !signe(g)? g: NULL;
229 5555830 : case t_INTMOD:
230 5555830 : return !signe(gel(g,2))? g: NULL;
231 3297 : case t_COMPLEX:
232 3297 : a = gisexactzero(gel(g,1)); if (!a) return NULL;
233 1281 : b = gisexactzero(gel(g,2)); if (!b) return NULL;
234 0 : return ggcd(a,b);
235 29309 : case t_FFELT:
236 29309 : return FF_equal0(g)? g: NULL;
237 581 : case t_QUAD:
238 581 : a = gisexactzero(gel(g,2)); if (!a) return NULL;
239 77 : b = gisexactzero(gel(g,3)); if (!b) return NULL;
240 7 : return ggcd(a,b);
241 15744 : case t_POLMOD:
242 15744 : return gisexactzero(gel(g,2));
243 25495611 : case t_POL:
244 25495611 : lx = lg(g); /* cater for Mod(0,2)*x^0 */
245 25495611 : if (lx == 2) return gen_0;
246 20163604 : if (lx == 3) return gisexactzero(gel(g,2));
247 16480922 : return NULL;
248 1191287 : case t_RFRAC:
249 1191287 : return gisexactzero(gel(g,1)); /* may occur: Mod(0,2)/x */
250 0 : case t_VEC: case t_COL: case t_MAT:
251 0 : a = gen_0;
252 0 : for (i=lg(g)-1; i; i--)
253 : {
254 0 : b = gisexactzero(gel(g,i));
255 0 : if (!b) return NULL;
256 0 : a = ggcd(a, b);
257 : }
258 0 : return a;
259 : }
260 604568 : return NULL;
261 : }
262 :
263 : int
264 792793339 : isrationalzero(GEN g)
265 : {
266 : long i;
267 792793339 : switch (typ(g))
268 : {
269 453461667 : case t_INT:
270 453461667 : return !signe(g);
271 42414908 : case t_COMPLEX:
272 42414908 : return isintzero(gel(g,1)) && isintzero(gel(g,2));
273 1435 : case t_QUAD:
274 1435 : return isintzero(gel(g,2)) && isintzero(gel(g,3));
275 507946 : case t_POLMOD:
276 507946 : return isrationalzero(gel(g,2));
277 138349796 : case t_POL: return lg(g) == 2;
278 91 : case t_VEC: case t_COL: case t_MAT:
279 322 : for (i=lg(g)-1; i; i--)
280 231 : if (!isrationalzero(gel(g,i))) return 0;
281 91 : return 1;
282 : }
283 158057496 : return 0;
284 : }
285 :
286 : int
287 2541470532 : gequal0(GEN x)
288 : {
289 2541470532 : switch(typ(x))
290 : {
291 2370107525 : case t_INT: case t_REAL: case t_POL: case t_SER:
292 2370107525 : return !signe(x);
293 :
294 7096672 : case t_INTMOD:
295 7096672 : return !signe(gel(x,2));
296 :
297 3100800 : case t_FFELT:
298 3100800 : return FF_equal0(x);
299 :
300 109365196 : case t_COMPLEX:
301 : /* is 0 iff norm(x) would be 0 (can happen with Re(x) and Im(x) != 0
302 : * only if Re(x) and Im(x) are of type t_REAL). See mp.c:addrr().
303 : */
304 109365196 : if (gequal0(gel(x,1)))
305 : {
306 8663336 : if (gequal0(gel(x,2))) return 1;
307 8270337 : if (typ(gel(x,1))!=t_REAL || typ(gel(x,2))!=t_REAL) return 0;
308 387455 : return (expo(gel(x,1))>=expo(gel(x,2)));
309 : }
310 100701860 : if (gequal0(gel(x,2)))
311 : {
312 2642197 : if (typ(gel(x,1))!=t_REAL || typ(gel(x,2))!=t_REAL) return 0;
313 2558202 : return (expo(gel(x,2))>=expo(gel(x,1)));
314 : }
315 98059663 : return 0;
316 :
317 2046559 : case t_PADIC:
318 2046559 : return !signe(padic_u(x));
319 :
320 1729 : case t_QUAD:
321 1729 : return gequal0(gel(x,2)) && gequal0(gel(x,3));
322 :
323 8556429 : case t_POLMOD:
324 8556429 : return gequal0(gel(x,2));
325 :
326 6113966 : case t_RFRAC:
327 6113966 : return gequal0(gel(x,1));
328 :
329 9912734 : case t_VEC: case t_COL: case t_MAT:
330 : {
331 : long i;
332 23826730 : for (i=lg(x)-1; i; i--)
333 19975502 : if (!gequal0(gel(x,i))) return 0;
334 3851228 : return 1;
335 : }
336 : }
337 25168922 : return 0;
338 : }
339 :
340 : /* x a t_POL or t_SER, return 1 if test(coeff(X,d)) is true and
341 : * coeff(X,i) = 0 for all i != d. Return 0 (false) otherwise */
342 : static int
343 22631465 : is_monomial_test(GEN x, long d, int(*test)(GEN))
344 : {
345 22631465 : long i, l = lg(x);
346 22631465 : if (typ(x) == t_SER)
347 : { /* "0" * x^v * (1+O(x)) ? v <= 0 or null ring */
348 602 : if (l == 3 && isexactzero(gel(x,2))) return d >= 2 || test(gel(x,2));
349 553 : if (d < 2) return 0; /* v > 0 */
350 : }
351 22631234 : if (d >= l)
352 : {
353 2302410 : if (typ(x) == t_POL) return 0; /* l = 2 */
354 : /* t_SER, v = 2-d <= 0 */
355 56 : if (!signe(x)) return 1;
356 : }
357 20328824 : else if (!test(gel(x,d))) return 0;
358 7585364 : for (i = 2; i < l; i++) /* 2 <= d < l */
359 5174294 : if (i != d && !gequal0(gel(x,i))) return 0;
360 2411070 : return 1;
361 : }
362 : static int
363 294560 : col_test(GEN x, int(*test)(GEN))
364 : {
365 294560 : long i, l = lg(x);
366 294560 : if (l == 1 || !test(gel(x,1))) return 0;
367 13167 : for (i = 2; i < l; i++)
368 11186 : if (!gequal0(gel(x,i))) return 0;
369 1981 : return 1;
370 : }
371 : static int
372 16366 : mat_test(GEN x, int(*test)(GEN))
373 : {
374 16366 : long i, j, l = lg(x);
375 16366 : if (l == 1) return 1;
376 16352 : if (l != lgcols(x)) return 0;
377 52283 : for (i = 1; i < l; i++)
378 140973 : for (j = 1; j < l; j++)
379 105042 : if (i == j) {
380 36015 : if (!test(gcoeff(x,i,i))) return 0;
381 : } else {
382 69027 : if (!gequal0(gcoeff(x,i,j))) return 0;
383 : }
384 16268 : return 1;
385 : }
386 :
387 : /* returns 1 whenever x = 1, and 0 otherwise */
388 : int
389 342288599 : gequal1(GEN x)
390 : {
391 342288599 : switch(typ(x))
392 : {
393 315195177 : case t_INT:
394 315195177 : return equali1(x);
395 :
396 82950 : case t_REAL:
397 : {
398 82950 : long s = signe(x);
399 82950 : if (!s) return expo(x) >= 0;
400 82320 : return s > 0 ? absrnz_equal1(x): 0;
401 : }
402 345971 : case t_INTMOD:
403 345971 : return is_pm1(gel(x,2)) || is_pm1(gel(x,1));
404 382817 : case t_POLMOD:
405 382817 : return !degpol(gel(x,1)) || gequal1(gel(x,2));
406 :
407 17108 : case t_FFELT:
408 17108 : return FF_equal1(x);
409 :
410 1708471 : case t_FRAC:
411 1708471 : return 0;
412 :
413 25681 : case t_COMPLEX:
414 25681 : return gequal1(gel(x,1)) && gequal0(gel(x,2));
415 :
416 166267 : case t_PADIC:
417 166267 : if (!signe(padic_u(x))) return valp(x) <= 0;
418 166225 : return valp(x) == 0 && gequal1(padic_u(x));
419 :
420 42 : case t_QUAD:
421 42 : return gequal1(gel(x,2)) && gequal0(gel(x,3));
422 :
423 22630793 : case t_POL: return is_monomial_test(x, 2, &gequal1);
424 476 : case t_SER: return is_monomial_test(x, 2 - valser(x), &gequal1);
425 :
426 1030454 : case t_RFRAC: return gequal(gel(x,1), gel(x,2));
427 294511 : case t_COL: return col_test(x, &gequal1);
428 16254 : case t_MAT: return mat_test(x, &gequal1);
429 : }
430 391627 : return 0;
431 : }
432 :
433 : /* returns 1 whenever the x = -1, 0 otherwise */
434 : int
435 74234552 : gequalm1(GEN x)
436 : {
437 : pari_sp av;
438 : GEN t;
439 :
440 74234552 : switch(typ(x))
441 : {
442 74224594 : case t_INT:
443 74224594 : return equalim1(x);
444 :
445 2016 : case t_REAL:
446 : {
447 2016 : long s = signe(x);
448 2016 : if (!s) return expo(x) >= 0;
449 1477 : return s < 0 ? absrnz_equal1(x): 0;
450 : }
451 4659 : case t_INTMOD:
452 4659 : av = avma; return gc_bool(av, equalii(addui(1,gel(x,2)), gel(x,1)));
453 :
454 154 : case t_FRAC:
455 154 : return 0;
456 :
457 42 : case t_FFELT:
458 42 : return FF_equalm1(x);
459 :
460 2576 : case t_COMPLEX:
461 2576 : return gequalm1(gel(x,1)) && gequal0(gel(x,2));
462 :
463 7 : case t_QUAD:
464 7 : return gequalm1(gel(x,2)) && gequal0(gel(x,3));
465 :
466 49 : case t_PADIC:
467 49 : t = padic_u(x); if (!signe(t)) return valp(x) <= 0;
468 21 : av = avma; return gc_bool(av, !valp(x) && equalii(addui(1,t), gel(x,3)));
469 :
470 56 : case t_POLMOD:
471 56 : return !degpol(gel(x,1)) || gequalm1(gel(x,2));
472 :
473 70 : case t_POL: return is_monomial_test(x, 2, &gequalm1);
474 126 : case t_SER: return is_monomial_test(x, 2 - valser(x), &gequalm1);
475 :
476 28 : case t_RFRAC:
477 28 : av = avma; return gc_bool(av, gmequal_try(gel(x,1), gel(x,2)));
478 49 : case t_COL: return col_test(x, &gequalm1);
479 112 : case t_MAT: return mat_test(x, &gequalm1);
480 : }
481 14 : return 0;
482 : }
483 :
484 : int
485 1494619 : gequalX(GEN x) { return typ(x) == t_POL && lg(x) == 4
486 10004131 : && isintzero(gel(x,2)) && isint1(gel(x,3)); }
487 :
488 : static int
489 3681020 : cmp_str(const char *x, const char *y)
490 : {
491 3681020 : int f = strcmp(x, y);
492 : return f > 0? 1
493 3681020 : : f? -1: 0;
494 : }
495 :
496 : static int
497 39285798 : cmp_universal_rec(GEN x, GEN y, long i0)
498 : {
499 39285798 : long i, lx = lg(x), ly = lg(y);
500 39285798 : if (lx < ly) return -1;
501 39282834 : if (lx > ly) return 1;
502 69619759 : for (i = i0; i < lx; i++)
503 : {
504 60675036 : int f = cmp_universal(gel(x,i), gel(y,i));
505 60675036 : if (f) return f;
506 : }
507 8944723 : return 0;
508 : }
509 : /* Universal "meaningless" comparison function. Transitive, returns 0 iff
510 : * gidentical(x,y) */
511 : int
512 89012831 : cmp_universal(GEN x, GEN y)
513 : {
514 89012831 : long lx, ly, i, tx = typ(x), ty = typ(y);
515 :
516 89012831 : if (tx < ty) return -1;
517 88634407 : if (ty < tx) return 1;
518 88097636 : switch(tx)
519 : {
520 44194799 : case t_INT: return cmpii(x,y);
521 3680999 : case t_STR: return cmp_str(GSTR(x),GSTR(y));
522 935956 : case t_REAL:
523 : case t_VECSMALL:
524 935956 : lx = lg(x);
525 935956 : ly = lg(y);
526 935956 : if (lx < ly) return -1;
527 887684 : if (lx > ly) return 1;
528 3590341 : for (i = 1; i < lx; i++)
529 : {
530 3481904 : if (x[i] < y[i]) return -1;
531 3112823 : if (x[i] > y[i]) return 1;
532 : }
533 108437 : return 0;
534 :
535 771423 : case t_POL:
536 : {
537 771423 : long X = x[1] & (VARNBITS|SIGNBITS);
538 771423 : long Y = y[1] & (VARNBITS|SIGNBITS);
539 771423 : if (X < Y) return -1;
540 771402 : if (X > Y) return 1;
541 771353 : return cmp_universal_rec(x, y, 2);
542 : }
543 881076 : case t_SER:
544 : case t_FFELT:
545 : case t_CLOSURE:
546 881076 : if (x[1] < y[1]) return -1;
547 881069 : if (x[1] > y[1]) return 1;
548 881062 : return cmp_universal_rec(x, y, 2);
549 :
550 35 : case t_LIST:
551 : {
552 35 : long tx = list_typ(x), ty = list_typ(y);
553 : GEN vx, vy;
554 : pari_sp av;
555 35 : if (tx < ty) return -1;
556 35 : if (tx > ty) return 1;
557 35 : vx = list_data(x);
558 35 : vy = list_data(y);
559 35 : if (!vx) return vy? -1: 0;
560 35 : if (!vy) return 1;
561 35 : av = avma;
562 35 : if (tx == t_LIST_MAP)
563 : {
564 14 : vx = maptomat_shallow(x);
565 14 : vy = maptomat_shallow(y);
566 : }
567 35 : return gc_int(av, cmp_universal_rec(vx, vy, 1));
568 : }
569 37633348 : default:
570 37633348 : return cmp_universal_rec(x, y, lontyp[tx]);
571 : }
572 : }
573 :
574 : static int
575 4595760 : cmpfrac(GEN x, GEN y)
576 : {
577 4595760 : pari_sp av = avma;
578 4595760 : GEN a = gel(x,1), b = gel(x,2);
579 4595760 : GEN c = gel(y,1), d = gel(y,2);
580 4595760 : return gc_bool(av, cmpii(mulii(a, d), mulii(b, c)));
581 : }
582 : static int
583 455849 : cmpifrac(GEN a, GEN y)
584 : {
585 455849 : pari_sp av = avma;
586 455849 : GEN c = gel(y,1), d = gel(y,2);
587 455849 : return gc_int(av, cmpii(mulii(a, d), c));
588 : }
589 : static int
590 81756 : cmprfrac(GEN a, GEN y)
591 : {
592 81756 : pari_sp av = avma;
593 81756 : GEN c = gel(y,1), d = gel(y,2);
594 81756 : return gc_int(av, cmpri(mulri(a, d), c));
595 : }
596 : static int
597 224 : cmpgen(GEN x, GEN y)
598 : {
599 224 : pari_sp av = avma;
600 224 : return gc_int(av, gsigne(gsub(x,y)));
601 : }
602 :
603 : /* returns the sign of x - y when it makes sense. 0 otherwise */
604 : int
605 231756946 : gcmp(GEN x, GEN y)
606 : {
607 231756946 : long tx = typ(x), ty = typ(y);
608 :
609 231756946 : if (tx == ty) /* generic case */
610 229544266 : switch(tx)
611 : {
612 110739283 : case t_INT: return cmpii(x, y);
613 114101372 : case t_REAL: return cmprr(x, y);
614 4595760 : case t_FRAC: return cmpfrac(x, y);
615 70 : case t_QUAD: return cmpgen(x, y);
616 21 : case t_STR: return cmp_str(GSTR(x), GSTR(y));
617 107760 : case t_INFINITY:
618 : {
619 107760 : long sx = inf_get_sign(x), sy = inf_get_sign(y);
620 107760 : if (sx < sy) return -1;
621 42 : if (sx > sy) return 1;
622 14 : return 0;
623 : }
624 : }
625 2212680 : if (ty == t_INFINITY) return -inf_get_sign(y);
626 1012985 : switch(tx)
627 : {
628 557283 : case t_INT:
629 557283 : switch(ty)
630 : {
631 244340 : case t_REAL: return cmpir(x, y);
632 312929 : case t_FRAC: return cmpifrac(x, y);
633 7 : case t_QUAD: return cmpgen(x, y);
634 : }
635 7 : break;
636 236514 : case t_REAL:
637 236514 : switch(ty)
638 : {
639 199198 : case t_INT: return cmpri(x, y);
640 37295 : case t_FRAC: return cmprfrac(x, y);
641 14 : case t_QUAD: return cmpgen(x, y);
642 : }
643 7 : break;
644 187395 : case t_FRAC:
645 187395 : switch(ty)
646 : {
647 142920 : case t_INT: return -cmpifrac(y, x);
648 44461 : case t_REAL: return -cmprfrac(y, x);
649 7 : case t_QUAD: return cmpgen(x, y);
650 : }
651 7 : break;
652 126 : case t_QUAD:
653 126 : return cmpgen(x, y);
654 31660 : case t_INFINITY: return inf_get_sign(x);
655 : }
656 28 : pari_err_TYPE2("comparison",x,y);
657 : return 0;/*LCOV_EXCL_LINE*/
658 : }
659 :
660 : int
661 762273 : gcmpsg(long s, GEN y)
662 : {
663 762273 : switch(typ(y))
664 : {
665 12453 : case t_INT: return cmpsi(s,y);
666 744626 : case t_REAL: return cmpsr(s,y);
667 5194 : case t_FRAC: {
668 5194 : pari_sp av = avma;
669 5194 : return gc_int(av, cmpii(mulsi(s,gel(y,2)), gel(y,1)));
670 : }
671 0 : case t_QUAD: {
672 0 : pari_sp av = avma;
673 0 : return gc_int(av, gsigne(gsubsg(s, y)));
674 : }
675 0 : case t_INFINITY: return -inf_get_sign(y);
676 : }
677 0 : pari_err_TYPE2("comparison",stoi(s),y);
678 : return 0; /* LCOV_EXCL_LINE */
679 : }
680 :
681 : static long
682 3257649 : roughtype(GEN x)
683 : {
684 3257649 : switch(typ(x))
685 : {
686 2114 : case t_MAT: return t_MAT;
687 742740 : case t_VEC: case t_COL: return t_VEC;
688 1634764 : case t_VECSMALL: return t_VECSMALL;
689 878031 : default: return t_INT;
690 : }
691 : }
692 :
693 : static int lexcmpsg(long x, GEN y);
694 : static int
695 42 : lexcmpgs(GEN x, long y) { return -lexcmpsg(y,x); }
696 : /* lexcmp(stoi(x),y), y t_VEC/t_COL/t_MAT */
697 : static int
698 21 : lexcmp_s_matvec(long x, GEN y)
699 : {
700 : int fl;
701 21 : if (lg(y)==1) return 1;
702 14 : fl = lexcmpsg(x,gel(y,1));
703 14 : if (fl) return fl;
704 7 : return -1;
705 : }
706 : /* x a scalar, y a t_VEC/t_COL/t_MAT */
707 : static int
708 357 : lexcmp_scal_matvec(GEN x, GEN y)
709 : {
710 : int fl;
711 357 : if (lg(y)==1) return 1;
712 357 : fl = lexcmp(x,gel(y,1));
713 357 : if (fl) return fl;
714 7 : return -1;
715 : }
716 : /* x a scalar, y a t_VECSMALL */
717 : static int
718 42 : lexcmp_scal_vecsmall(GEN x, GEN y)
719 : {
720 : int fl;
721 42 : if (lg(y)==1) return 1;
722 42 : fl = lexcmpgs(x, y[1]);
723 42 : if (fl) return fl;
724 0 : return -1;
725 : }
726 :
727 : /* tx = ty = t_MAT, or x and y are both vect_t */
728 : static int
729 372161 : lexcmp_similar(GEN x, GEN y)
730 : {
731 372161 : long i, lx = lg(x), ly = lg(y), l = minss(lx,ly);
732 457176 : for (i=1; i<l; i++)
733 : {
734 425536 : int fl = lexcmp(gel(x,i),gel(y,i));
735 425536 : if (fl) return fl;
736 : }
737 31640 : if (lx == ly) return 0;
738 35 : return (lx < ly)? -1 : 1;
739 : }
740 : /* x a t_VECSMALL, y a t_VEC/t_COL ~ lexcmp_similar */
741 : static int
742 154 : lexcmp_vecsmall_vec(GEN x, GEN y)
743 : {
744 154 : long i, lx = lg(x), ly = lg(y), l = minss(lx,ly);
745 343 : for (i=1; i<l; i++)
746 : {
747 287 : int fl = lexcmpsg(x[i], gel(y,i));
748 287 : if (fl) return fl;
749 : }
750 56 : if (lx == ly) return 0;
751 21 : return (lx < ly)? -1 : 1;
752 : }
753 :
754 : /* x t_VEC/t_COL, y t_MAT */
755 : static int
756 98 : lexcmp_vec_mat(GEN x, GEN y)
757 : {
758 : int fl;
759 98 : if (lg(x)==1) return -1;
760 98 : if (lg(y)==1) return 1;
761 98 : fl = lexcmp_similar(x,gel(y,1));
762 98 : if (fl) return fl;
763 7 : return -1;
764 : }
765 : /* x t_VECSMALl, y t_MAT ~ lexcmp_vec_mat */
766 : static int
767 42 : lexcmp_vecsmall_mat(GEN x, GEN y)
768 : {
769 : int fl;
770 42 : if (lg(x)==1) return -1;
771 42 : if (lg(y)==1) return 1;
772 42 : fl = lexcmp_vecsmall_vec(x, gel(y,1));
773 42 : if (fl) return fl;
774 0 : return -1;
775 : }
776 :
777 : /* x a t_VECSMALL, not y */
778 : static int
779 196 : lexcmp_vecsmall_other(GEN x, GEN y, long ty)
780 : {
781 196 : switch(ty)
782 : {
783 42 : case t_MAT: return lexcmp_vecsmall_mat(x, y);
784 112 : case t_VEC: return lexcmp_vecsmall_vec(x, y);
785 42 : default: return -lexcmp_scal_vecsmall(y, x); /*y scalar*/
786 : }
787 : }
788 :
789 : /* lexcmp(stoi(s), y) */
790 : static int
791 343 : lexcmpsg(long x, GEN y)
792 : {
793 343 : switch(roughtype(y))
794 : {
795 21 : case t_MAT:
796 : case t_VEC:
797 21 : return lexcmp_s_matvec(x,y);
798 14 : case t_VECSMALL: /* ~ lexcmp_scal_matvec */
799 14 : if (lg(y)==1) return 1;
800 7 : return (x > y[1])? 1: -1;
801 308 : default: return gcmpsg(x,y);
802 : }
803 : }
804 :
805 : /* as gcmp for vector/matrices, using lexicographic ordering on components */
806 : static int
807 1628653 : lexcmp_i(GEN x, GEN y)
808 : {
809 1628653 : const long tx = roughtype(x), ty = roughtype(y);
810 1628653 : if (tx == ty)
811 1628002 : switch(tx)
812 : {
813 372063 : case t_MAT:
814 372063 : case t_VEC: return lexcmp_similar(x,y);
815 817277 : case t_VECSMALL: return vecsmall_lexcmp(x,y);
816 438662 : default: return gcmp(x,y);
817 : }
818 651 : if (tx == t_VECSMALL) return lexcmp_vecsmall_other(x,y,ty);
819 518 : if (ty == t_VECSMALL) return -lexcmp_vecsmall_other(y,x,tx);
820 :
821 455 : if (tx == t_INT) return lexcmp_scal_matvec(x,y); /*scalar*/
822 203 : if (ty == t_INT) return -lexcmp_scal_matvec(y,x);
823 :
824 98 : if (ty==t_MAT) return lexcmp_vec_mat(x,y);
825 42 : return -lexcmp_vec_mat(y,x); /*tx==t_MAT*/
826 : }
827 : int
828 1628653 : lexcmp(GEN x, GEN y)
829 : {
830 1628653 : pari_sp av = avma;
831 1628653 : if (typ(x) == t_COMPLEX)
832 : {
833 1001 : x = mkvec2(gel(x,1), gel(x,2));
834 1001 : if (typ(y) == t_COMPLEX) y = mkvec2(gel(y,1), gel(y,2));
835 49 : else y = mkvec2(y, gen_0);
836 : }
837 1627652 : else if (typ(y) == t_COMPLEX)
838 : {
839 63 : x = mkvec2(x, gen_0);
840 63 : y = mkvec2(gel(y,1), gel(y,2));
841 : }
842 1628653 : return gc_int(av, lexcmp_i(x, y));
843 : }
844 :
845 : /*****************************************************************/
846 : /* */
847 : /* EQUALITY */
848 : /* returns 1 if x == y, 0 otherwise */
849 : /* */
850 : /*****************************************************************/
851 : /* x,y t_POL */
852 : static int
853 5634065 : polidentical(GEN x, GEN y)
854 : {
855 : long lx;
856 5634065 : if (x[1] != y[1]) return 0;
857 5633967 : lx = lg(x); if (lg(y) != lg(x)) return 0;
858 37916365 : for (lx--; lx >= 2; lx--) if (!gidentical(gel(x,lx), gel(y,lx))) return 0;
859 5633470 : return 1;
860 : }
861 : /* x,y t_SER */
862 : static int
863 14 : seridentical(GEN x, GEN y) { return polidentical(x,y); }
864 : /* typ(x) = typ(y) = t_VEC/COL/MAT */
865 : static int
866 5440915 : vecidentical(GEN x, GEN y)
867 : {
868 : long i;
869 5440915 : if ((x[0] ^ y[0]) & (TYPBITS|LGBITS)) return 0;
870 17253355 : for (i = lg(x)-1; i; i--)
871 13123596 : if (! gidentical(gel(x,i),gel(y,i)) ) return 0;
872 4129759 : return 1;
873 : }
874 : static int
875 1547 : identicalrr(GEN x, GEN y)
876 : {
877 1547 : long i, lx = lg(x);
878 1547 : if (lg(y) != lx) return 0;
879 1547 : if (x[1] != y[1]) return 0;
880 1540 : if (!signe(x)) return 1;
881 4924 : i=2; while (i<lx && x[i]==y[i]) i++;
882 1407 : return (i == lx);
883 : }
884 :
885 : static int
886 70 : closure_identical(GEN x, GEN y)
887 : {
888 70 : if (lg(x)!=lg(y) || x[1]!=y[1]) return 0;
889 56 : if (!gidentical(gel(x,2),gel(y,2)) || !gidentical(gel(x,3),gel(y,3))
890 56 : || !gidentical(gel(x,4),gel(y,4))) return 0;
891 42 : if (lg(x)<8) return 1;
892 0 : return gidentical(gel(x,7),gel(y,7));
893 : }
894 :
895 : static int
896 343 : list_cmp(GEN x, GEN y, int cmp(GEN x, GEN y))
897 : {
898 343 : int t = list_typ(x);
899 : GEN vx, vy;
900 : long lvx, lvy;
901 343 : if (list_typ(y)!=t) return 0;
902 343 : vx = list_data(x);
903 343 : vy = list_data(y);
904 343 : lvx = vx ? lg(vx): 1;
905 343 : lvy = vy ? lg(vy): 1;
906 343 : if (lvx==1 && lvy==1) return 1;
907 329 : if (lvx != lvy) return 0;
908 301 : switch (t)
909 : {
910 280 : case t_LIST_MAP:
911 : {
912 280 : pari_sp av = avma;
913 280 : GEN mx = maptomat_shallow(x), my = maptomat_shallow(y);
914 280 : int ret = gidentical(gel(mx, 1), gel(my, 1)) && cmp(gel(mx, 2), gel(my, 2));
915 280 : return gc_bool(av, ret);
916 : }
917 21 : default:
918 21 : return cmp(vx, vy);
919 : }
920 : }
921 :
922 : int
923 123903204 : gidentical(GEN x, GEN y)
924 : {
925 : long tx;
926 :
927 123903204 : if (x == y) return 1;
928 118396885 : tx = typ(x); if (typ(y) != tx) return 0;
929 118161241 : switch(tx)
930 : {
931 39946240 : case t_INT:
932 39946240 : return equalii(x,y);
933 :
934 1547 : case t_REAL:
935 1547 : return identicalrr(x,y);
936 :
937 426602 : case t_FRAC: case t_INTMOD:
938 426602 : return equalii(gel(x,2), gel(y,2)) && equalii(gel(x,1), gel(y,1));
939 :
940 343 : case t_COMPLEX:
941 343 : return gidentical(gel(x,2),gel(y,2)) && gidentical(gel(x,1),gel(y,1));
942 21 : case t_PADIC:
943 21 : return valp(x) == valp(y) && precp(x) == precp(y)
944 14 : && equalii(padic_p(x), padic_p(y))
945 42 : && equalii(padic_u(x), padic_u(y));
946 3892 : case t_POLMOD:
947 3892 : return gidentical(gel(x,2),gel(y,2)) && polidentical(gel(x,1),gel(y,1));
948 5634030 : case t_POL:
949 5634030 : return polidentical(x,y);
950 14 : case t_SER:
951 14 : return seridentical(x,y);
952 5320 : case t_FFELT:
953 5320 : return FF_equal(x,y);
954 :
955 423760 : case t_QFB:
956 423760 : return equalii(gel(x,1),gel(y,1))
957 423753 : && equalii(gel(x,2),gel(y,2))
958 847513 : && equalii(gel(x,3),gel(y,3));
959 :
960 14 : case t_QUAD:
961 14 : return ZX_equal(gel(x,1),gel(y,1))
962 7 : && gidentical(gel(x,2),gel(y,2))
963 21 : && gidentical(gel(x,3),gel(y,3));
964 :
965 7 : case t_RFRAC:
966 7 : return gidentical(gel(x,1),gel(y,1)) && gidentical(gel(x,2),gel(y,2));
967 :
968 70 : case t_STR:
969 70 : return !strcmp(GSTR(x),GSTR(y));
970 5440915 : case t_VEC: case t_COL: case t_MAT:
971 5440915 : return vecidentical(x,y);
972 66278256 : case t_VECSMALL:
973 66278256 : return zv_equal(x,y);
974 28 : case t_CLOSURE:
975 28 : return closure_identical(x,y);
976 161 : case t_LIST:
977 161 : return list_cmp(x, y, gidentical);
978 21 : case t_INFINITY: return gidentical(gel(x,1),gel(y,1));
979 : }
980 0 : return 0;
981 : }
982 : /* x,y t_POL in the same variable */
983 : static int
984 12918186 : polequal(GEN x, GEN y)
985 : {
986 : long lx, ly;
987 : /* Can't do that: Mod(0,1)*x^0 == x^0
988 : if (signe(x) != signe(y)) return 0; */
989 12918186 : lx = lg(x); ly = lg(y);
990 12918186 : while (lx > ly) if (!gequal0(gel(x,--lx))) return 0;
991 12915701 : while (ly > lx) if (!gequal0(gel(y,--ly))) return 0;
992 45045121 : for (lx--; lx >= 2; lx--) if (!gequal(gel(x,lx), gel(y,lx))) return 0;
993 12844418 : return 1;
994 : }
995 :
996 : /* x,y t_SER in the same variable */
997 : static int
998 420 : serequal(GEN x, GEN y)
999 : {
1000 : long LX, LY, lx, ly, vx, vy;
1001 420 : if (!signe(x) && !signe(y)) return 1;
1002 56 : lx = lg(x); vx = valser(x); LX = lx + vx;
1003 56 : ly = lg(y); vy = valser(y); LY = ly + vy;
1004 56 : if (LX > LY) lx = LY - vx; else ly = LX - vy;
1005 282877 : while (lx >= 3 && ly >= 3)
1006 282821 : if (!gequal(gel(x,--lx), gel(y,--ly))) return 0;
1007 56 : while(--ly >= 2) if (!gequal0(gel(y,ly))) return 0;
1008 84 : while(--lx >= 2) if (!gequal0(gel(x,lx))) return 0;
1009 49 : return 1;
1010 : }
1011 :
1012 : /* typ(x) = typ(y) = t_VEC/COL/MAT */
1013 : static int
1014 5577413 : vecequal(GEN x, GEN y)
1015 : {
1016 : long i;
1017 5577413 : if ((x[0] ^ y[0]) & (TYPBITS|LGBITS)) return 0;
1018 18606609 : for (i = lg(x)-1; i; i--)
1019 16242572 : if (! gequal(gel(x,i),gel(y,i)) ) return 0;
1020 2364037 : return 1;
1021 : }
1022 :
1023 : int
1024 244360168 : gequal(GEN x, GEN y)
1025 : {
1026 : pari_sp av;
1027 : GEN A, B, a, b;
1028 : long tx, ty;
1029 :
1030 244360168 : if (x == y) return 1;
1031 210321187 : tx = typ(x); ty = typ(y);
1032 210321187 : if (tx == ty)
1033 202489556 : switch(tx)
1034 : {
1035 172055844 : case t_INT:
1036 172055844 : return equalii(x,y);
1037 :
1038 20488 : case t_REAL:
1039 20488 : return equalrr(x,y);
1040 :
1041 299284 : case t_FRAC:
1042 299284 : return equalii(gel(x,2), gel(y,2)) && equalii(gel(x,1), gel(y,1));
1043 :
1044 6510007 : case t_INTMOD:
1045 6510007 : A = gel(x,1); B = gel(y,1);
1046 6510007 : a = gel(x,2); b = gel(y,2);
1047 6510007 : if (equalii(A, B)) return equalii(a, b);
1048 14 : av = avma; A = gcdii(A, B);
1049 14 : return gc_bool(av, equalii(modii(a,A), modii(b,A)));
1050 :
1051 1316 : case t_COMPLEX:
1052 1316 : return gequal(gel(x,2),gel(y,2)) && gequal(gel(x,1),gel(y,1));
1053 770 : case t_PADIC:
1054 770 : if (!equalii(padic_p(x), padic_p(y))) return 0;
1055 770 : av = avma; return gc_bool(av, gequal0(gsub(x,y)));
1056 :
1057 3189703 : case t_POLMOD:
1058 3189703 : A = gel(x,1); B = gel(y,1);
1059 3189703 : if (varn(A) != varn(B)) break;
1060 3189682 : a = gel(x,2); b = gel(y,2);
1061 3189682 : if (RgX_equal_var(A, B)) return gequal(a,b);
1062 14 : av = avma; A = ggcd(A, B);
1063 14 : return gc_bool(av, gequal(gmod(a,A), gmod(b,A)));
1064 :
1065 12930338 : case t_POL:
1066 12930338 : if (varn(x) != varn(y)) break;
1067 12918186 : return polequal(x,y);
1068 420 : case t_SER:
1069 420 : if (varn(x) != varn(y)) break;
1070 420 : return serequal(x,y);
1071 :
1072 250713 : case t_FFELT:
1073 250713 : return FF_equal(x,y);
1074 :
1075 91513 : case t_QFB:
1076 91513 : return equalii(gel(x,1),gel(y,1))
1077 72986 : && equalii(gel(x,2),gel(y,2))
1078 164499 : && equalii(gel(x,3),gel(y,3));
1079 :
1080 84 : case t_QUAD:
1081 84 : return ZX_equal(gel(x,1),gel(y,1))
1082 77 : && gequal(gel(x,2),gel(y,2))
1083 161 : && gequal(gel(x,3),gel(y,3));
1084 :
1085 73759 : case t_RFRAC:
1086 : {
1087 73759 : GEN a = gel(x,1), b = gel(x,2), c = gel(y,1), d = gel(y,2);
1088 73759 : if (gequal(b,d)) return gequal(a,c); /* simple case */
1089 0 : av = avma;
1090 0 : a = simplify_shallow(gmul(a,d));
1091 0 : b = simplify_shallow(gmul(b,c));
1092 0 : return gc_bool(av, gequal(a,b));
1093 : }
1094 :
1095 65808 : case t_STR:
1096 65808 : return !strcmp(GSTR(x),GSTR(y));
1097 5577413 : case t_VEC: case t_COL: case t_MAT:
1098 5577413 : return vecequal(x,y);
1099 1421844 : case t_VECSMALL:
1100 1421844 : return zv_equal(x,y);
1101 182 : case t_LIST:
1102 182 : return list_cmp(x, y, gequal);
1103 42 : case t_CLOSURE:
1104 42 : return closure_identical(x,y);
1105 28 : case t_INFINITY:
1106 28 : return gequal(gel(x,1),gel(y,1));
1107 : }
1108 7843804 : if (is_noncalc_t(tx) || is_noncalc_t(ty)) return 0;
1109 7843720 : if (tx == t_INT && !signe(x)) return gequal0(y);
1110 7800869 : if (ty == t_INT && !signe(y)) return gequal0(x);
1111 3195604 : (void)&av; av = avma; /* emulate volatile */
1112 3195604 : return gc_bool(av, gequal_try(x, y));
1113 : }
1114 :
1115 : int
1116 44058 : gequalsg(long s, GEN x)
1117 44058 : { pari_sp av = avma; return gc_bool(av, gequal(stoi(s), x)); }
1118 :
1119 : /* a and b are t_INT, t_FRAC, t_REAL or t_COMPLEX of those. Check whether
1120 : * a-b is invertible */
1121 : int
1122 2499437 : cx_approx_equal(GEN a, GEN b)
1123 : {
1124 2499437 : pari_sp av = avma;
1125 : GEN d;
1126 2499437 : if (a == b) return 1;
1127 639128 : d = gsub(a,b);
1128 639128 : return gc_bool(av, gequal0(d) || (typ(d)==t_COMPLEX && gequal0(cxnorm(d))));
1129 : }
1130 : static int
1131 1808884 : r_approx0(GEN x, long e) { return e - expo(x) > bit_prec(x); }
1132 : /* x ~ 0 compared to reference y */
1133 : int
1134 2506963 : cx_approx0(GEN x, GEN y)
1135 : {
1136 : GEN a, b;
1137 : long e;
1138 2506963 : switch(typ(x))
1139 : {
1140 511 : case t_COMPLEX:
1141 511 : a = gel(x,1); b = gel(x,2);
1142 511 : if (typ(a) != t_REAL)
1143 : {
1144 14 : if (!gequal0(a)) return 0;
1145 0 : a = NULL;
1146 : }
1147 497 : else if (!signe(a)) a = NULL;
1148 497 : if (typ(b) != t_REAL)
1149 : {
1150 0 : if (!gequal0(b)) return 0;
1151 0 : if (!a) return 1;
1152 0 : b = NULL;
1153 : }
1154 497 : else if (!signe(b))
1155 : {
1156 28 : if (!a) return 1;
1157 7 : b = NULL;
1158 : }
1159 : /* a or b is != NULL iff it is non-zero t_REAL; one of them is */
1160 476 : e = gexpo(y);
1161 476 : return (!a || r_approx0(a, e)) && (!b || r_approx0(b, e));
1162 1808597 : case t_REAL:
1163 1808597 : return !signe(x) || r_approx0(x, gexpo(y));
1164 697855 : default:
1165 697855 : return gequal0(x);
1166 : }
1167 : }
1168 : /*******************************************************************/
1169 : /* */
1170 : /* VALUATION */
1171 : /* p is either a t_INT or a t_POL. */
1172 : /* returns the largest exponent of p dividing x when this makes */
1173 : /* sense : error for types real, integermod and polymod if p does */
1174 : /* not divide the modulus, q-adic if q!=p. */
1175 : /* */
1176 : /*******************************************************************/
1177 :
1178 : static long
1179 137186 : minval(GEN x, GEN p)
1180 : {
1181 137186 : long i,k, val = LONG_MAX, lx = lg(x);
1182 376943 : for (i=lontyp[typ(x)]; i<lx; i++)
1183 : {
1184 239757 : k = gvaluation(gel(x,i),p);
1185 239757 : if (k < val) val = k;
1186 : }
1187 137186 : return val;
1188 : }
1189 :
1190 : static int
1191 91 : intdvd(GEN x, GEN y, GEN *z) { GEN r; *z = dvmdii(x,y,&r); return (r==gen_0); }
1192 :
1193 : /* x t_FRAC, p t_INT, return v_p(x) */
1194 : static long
1195 292596 : frac_val(GEN x, GEN p) {
1196 292596 : long v = Z_pval(gel(x,2),p);
1197 292596 : if (v) return -v;
1198 292454 : return Z_pval(gel(x,1),p);
1199 : }
1200 : long
1201 9427685 : Q_pval(GEN x, GEN p)
1202 : {
1203 9427685 : if (lgefint(p) == 3) return Q_lval(x, uel(p,2));
1204 555 : return (typ(x)==t_INT)? Z_pval(x, p): frac_val(x, p);
1205 : }
1206 :
1207 : static long
1208 381841 : frac_lval(GEN x, ulong p) {
1209 381841 : long v = Z_lval(gel(x,2),p);
1210 381841 : if (v) return -v;
1211 232161 : return Z_lval(gel(x,1),p);
1212 : }
1213 : long
1214 9431876 : Q_lval(GEN x, ulong p){return (typ(x)==t_INT)? Z_lval(x, p): frac_lval(x, p);}
1215 :
1216 : long
1217 6304619 : Q_pvalrem(GEN x, GEN p, GEN *y)
1218 : {
1219 : GEN a, b;
1220 : long v;
1221 6304619 : if (lgefint(p) == 3) return Q_lvalrem(x, uel(p,2), y);
1222 5810 : if (typ(x) == t_INT) return Z_pvalrem(x, p, y);
1223 0 : a = gel(x,1);
1224 0 : b = gel(x,2);
1225 0 : v = Z_pvalrem(b, p, &b);
1226 0 : if (v) { *y = isint1(b)? a: mkfrac(a, b); return -v; }
1227 0 : v = Z_pvalrem(a, p, &a);
1228 0 : *y = mkfrac(a, b); return v;
1229 : }
1230 : long
1231 6303023 : Q_lvalrem(GEN x, ulong p, GEN *y)
1232 : {
1233 : GEN a, b;
1234 : long v;
1235 6303023 : if (typ(x) == t_INT) return Z_lvalrem(x, p, y);
1236 398420 : a = gel(x,1);
1237 398420 : b = gel(x,2);
1238 398420 : v = Z_lvalrem(b, p, &b);
1239 398420 : if (v) { *y = isint1(b)? a: mkfrac(a, b); return -v; }
1240 228917 : v = Z_lvalrem(a, p, &a);
1241 228917 : *y = mkfrac(a, b); return v;
1242 : }
1243 :
1244 : long
1245 1401625 : gvaluation(GEN x, GEN p)
1246 : {
1247 1401625 : long tx = typ(x), tp;
1248 : pari_sp av;
1249 :
1250 1401625 : if (!p)
1251 28 : switch(tx)
1252 : {
1253 7 : case t_PADIC: return valp(x);
1254 7 : case t_POL: return RgX_val(x);
1255 7 : case t_SER: return valser(x);
1256 7 : default: pari_err_TYPE("gvaluation", x);
1257 : }
1258 1401597 : tp = typ(p);
1259 1401597 : switch(tp)
1260 : {
1261 1341698 : case t_INT:
1262 1341698 : if (signe(p) && !is_pm1(p)) break;
1263 28 : pari_err_DOMAIN("gvaluation", "p", "=", p, p);
1264 59892 : case t_POL:
1265 59892 : if (degpol(p) > 0) break;
1266 : default:
1267 7 : pari_err_DOMAIN("gvaluation", "p", "=", p, p);
1268 : }
1269 :
1270 1401562 : switch(tx)
1271 : {
1272 224941 : case t_INT:
1273 224941 : if (!signe(x)) return LONG_MAX;
1274 160758 : if (tp == t_POL) return 0;
1275 160142 : return Z_pval(x,p);
1276 :
1277 49 : case t_REAL:
1278 49 : if (tp == t_POL) return 0;
1279 21 : break;
1280 :
1281 28 : case t_FFELT:
1282 28 : if (tp == t_POL) return FF_equal0(x)? LONG_MAX: 0;
1283 14 : break;
1284 :
1285 105 : case t_INTMOD: {
1286 105 : GEN a = gel(x,1), b = gel(x,2);
1287 : long val;
1288 133 : if (tp == t_POL) return signe(b)? 0: LONG_MAX;
1289 42 : av = avma;
1290 42 : if (!intdvd(a, p, &a)) break;
1291 28 : if (!intdvd(b, p, &b)) return gc_long(av,0);
1292 14 : val = 1; while (intdvd(a,p,&a) && intdvd(b,p,&b)) val++;
1293 14 : return gc_long(av,val);
1294 : }
1295 :
1296 292516 : case t_FRAC:
1297 292516 : if (tp == t_POL) return 0;
1298 292502 : return frac_val(x, p);
1299 :
1300 721908 : case t_PADIC:
1301 721908 : if (tp == t_POL) return 0;
1302 721887 : if (!equalii(p, padic_p(x))) break;
1303 721880 : return valp(x);
1304 :
1305 35 : case t_POLMOD: {
1306 35 : GEN a = gel(x,1), b = gel(x,2);
1307 : long v, val;
1308 35 : if (tp == t_INT) return gvaluation(b,p);
1309 21 : v = varn(p);
1310 21 : if (varn(a) != v) return 0;
1311 21 : av = avma;
1312 21 : a = RgX_divrem(a, p, ONLY_DIVIDES);
1313 21 : if (!a) break;
1314 28 : if (typ(b) != t_POL || varn(b) != v ||
1315 21 : !(b = RgX_divrem(b, p, ONLY_DIVIDES)) ) return gc_long(av,0);
1316 7 : val = 1;
1317 28 : while ((a = RgX_divrem(a, p, ONLY_DIVIDES)) &&
1318 21 : (b = RgX_divrem(b, p, ONLY_DIVIDES)) ) val++;
1319 7 : return gc_long(av,val);
1320 : }
1321 161315 : case t_POL: {
1322 161315 : if (tp == t_POL) {
1323 58023 : long vp = varn(p), vx = varn(x);
1324 58023 : if (vp == vx)
1325 : {
1326 : long val;
1327 24262 : if (RgX_is_monomial(p))
1328 : {
1329 24227 : val = RgX_val(x); if (val == LONG_MAX) return LONG_MAX;
1330 13384 : return val / degpol(p);
1331 : }
1332 35 : if (!signe(x)) return LONG_MAX;
1333 21 : av = avma;
1334 21 : for (val=0; ; val++)
1335 : {
1336 35 : x = RgX_divrem(x,p,ONLY_DIVIDES);
1337 35 : if (!x) return gc_long(av,val);
1338 14 : if (gc_needed(av,1))
1339 : {
1340 0 : if(DEBUGMEM>1) pari_warn(warnmem,"gvaluation");
1341 0 : x = gc_GEN(av, x);
1342 : }
1343 : }
1344 : }
1345 33761 : if (varncmp(vx, vp) > 0) return 0;
1346 : }
1347 137046 : return minval(x,p);
1348 : }
1349 :
1350 490 : case t_SER: {
1351 490 : if (tp == t_POL) {
1352 476 : long vp = varn(p), vx = varn(x);
1353 476 : if (vp == vx)
1354 : {
1355 469 : long val = RgX_val(p);
1356 469 : if (!val) pari_err_DOMAIN("gvaluation", "p", "=", p, p);
1357 462 : return (long)(valser(x) / val);
1358 : }
1359 7 : if (varncmp(vx, vp) > 0) return 0;
1360 : }
1361 14 : return minval(x,p);
1362 : }
1363 :
1364 49 : case t_RFRAC:
1365 49 : return gvaluation(gel(x,1),p) - gvaluation(gel(x,2),p);
1366 :
1367 126 : case t_COMPLEX: case t_QUAD: case t_VEC: case t_COL: case t_MAT:
1368 126 : return minval(x,p);
1369 : }
1370 63 : pari_err_OP("valuation", x,p);
1371 : return 0; /* LCOV_EXCL_LINE */
1372 : }
1373 : GEN
1374 3934 : gpvaluation(GEN x, GEN p)
1375 : {
1376 3934 : long v = gvaluation(x,p);
1377 3829 : return v == LONG_MAX? mkoo(): stoi(v);
1378 : }
1379 :
1380 : /* x is nonzero */
1381 : long
1382 90418412 : u_lvalrem(ulong x, ulong p, ulong *py)
1383 : {
1384 : ulong vx;
1385 90418412 : if (p == 2) { vx = vals(x); *py = x >> vx; return vx; }
1386 80153289 : for(vx = 0;;)
1387 : {
1388 131122485 : if (x % p) { *py = x; return vx; }
1389 50969196 : x /= p; /* gcc is smart enough to make a single div */
1390 50969196 : vx++;
1391 : }
1392 : }
1393 : long
1394 71118590 : u_lval(ulong x, ulong p)
1395 : {
1396 : ulong vx;
1397 71118590 : if (p == 2) return vals(x);
1398 67729141 : for(vx = 0;;)
1399 : {
1400 115948724 : if (x % p) return vx;
1401 48219583 : x /= p; /* gcc is smart enough to make a single div */
1402 48219583 : vx++;
1403 : }
1404 : }
1405 :
1406 : long
1407 1839888 : z_lval(long s, ulong p) { return u_lval(labs(s), p); }
1408 : long
1409 87353 : z_lvalrem(long s, ulong p, long *py)
1410 : {
1411 : long v;
1412 87353 : if (s < 0)
1413 : {
1414 0 : ulong u = (ulong)-s;
1415 0 : v = u_lvalrem(u, p, &u);
1416 0 : *py = -(long)u;
1417 : }
1418 : else
1419 : {
1420 87353 : ulong u = (ulong)s;
1421 87353 : v = u_lvalrem(u, p, &u);
1422 87353 : *py = (long)u;
1423 : }
1424 87353 : return v;
1425 : }
1426 : /* assume |p| > 1 */
1427 : long
1428 1318037 : z_pval(long s, GEN p)
1429 : {
1430 1318037 : if (lgefint(p) > 3) return 0;
1431 1318037 : return z_lval(s, uel(p,2));
1432 : }
1433 : /* assume |p| > 1 */
1434 : long
1435 399 : z_pvalrem(long s, GEN p, long *py)
1436 : {
1437 399 : if (lgefint(p) > 3) { *py = s; return 0; }
1438 399 : return z_lvalrem(s, uel(p,2), py);
1439 : }
1440 :
1441 : /* return v_q(x) and set *py = x / q^v_q(x), using divide & conquer */
1442 : static long
1443 2159121 : Z_pvalrem_DC(GEN x, GEN q, GEN *py)
1444 : {
1445 2159121 : GEN r, z = dvmdii(x, q, &r);
1446 : long v;
1447 2159121 : if (r != gen_0) { *py = x; return 0; }
1448 1488773 : if (2 * lgefint(q) <= lgefint(z)+3) /* avoid squaring if pointless */
1449 1472625 : v = Z_pvalrem_DC(z, sqri(q), py) << 1;
1450 : else
1451 16148 : { v = 0; *py = z; }
1452 1488773 : z = dvmdii(*py, q, &r);
1453 1488773 : if (r != gen_0) return v + 1;
1454 622531 : *py = z; return v + 2;
1455 : }
1456 :
1457 : static const long VAL_DC_THRESHOLD = 16;
1458 :
1459 : long
1460 68998038 : Z_lval(GEN x, ulong p)
1461 : {
1462 : long vx;
1463 : pari_sp av;
1464 68998038 : if (p == 2) return vali(x);
1465 53758182 : if (lgefint(x) == 3) return u_lval(uel(x,2), p);
1466 2758087 : av = avma;
1467 2758087 : for(vx = 0;;)
1468 11887194 : {
1469 : ulong r;
1470 14645281 : GEN q = absdiviu_rem(x, p, &r);
1471 14645281 : if (r) break;
1472 12075973 : vx++; x = q;
1473 12075973 : if (vx == VAL_DC_THRESHOLD) {
1474 188779 : if (p == 1) pari_err_DOMAIN("Z_lval", "p", "=", gen_1, gen_1);
1475 188779 : vx += Z_pvalrem_DC(x, sqru(p), &x) << 1;
1476 188779 : q = absdiviu_rem(x, p, &r); if (!r) vx++;
1477 188779 : break;
1478 : }
1479 : }
1480 2758087 : return gc_long(av,vx);
1481 : }
1482 : long
1483 65075951 : Z_lvalrem(GEN x, ulong p, GEN *py)
1484 : {
1485 : long vx, sx;
1486 : pari_sp av;
1487 65075951 : if (p == 2) { vx = vali(x); *py = shifti(x, -vx); return vx; }
1488 51176975 : if (lgefint(x) == 3) {
1489 : ulong u;
1490 44405194 : vx = u_lvalrem(uel(x,2), p, &u);
1491 44405194 : *py = signe(x) < 0? utoineg(u): utoipos(u);
1492 44405194 : return vx;
1493 : }
1494 6771781 : av = avma; (void)new_chunk(lgefint(x));
1495 6771781 : sx = signe(x);
1496 6771781 : for(vx = 0;;)
1497 19141282 : {
1498 : ulong r;
1499 25913063 : GEN q = absdiviu_rem(x, p, &r);
1500 25913063 : if (r) break;
1501 19638370 : vx++; x = q;
1502 19638370 : if (vx == VAL_DC_THRESHOLD) {
1503 497088 : if (p == 1) pari_err_DOMAIN("Z_lvalrem", "p", "=", gen_1, gen_1);
1504 497088 : vx += Z_pvalrem_DC(x, sqru(p), &x) << 1;
1505 497088 : q = absdiviu_rem(x, p, &r); if (!r) { vx++; x = q; }
1506 497088 : break;
1507 : }
1508 : }
1509 6771781 : set_avma(av); *py = icopy(x); setsigne(*py, sx); return vx;
1510 : }
1511 :
1512 : /* Is |q| <= p ? */
1513 : static int
1514 15096244 : isless_iu(GEN q, ulong p) {
1515 15096244 : long l = lgefint(q);
1516 15096244 : return l==2 || (l == 3 && uel(q,2) <= p);
1517 : }
1518 :
1519 : long
1520 136147098 : u_lvalrem_stop(ulong *n, ulong p, int *stop)
1521 : {
1522 136147098 : ulong N = *n, q = N / p, r = N % p; /* gcc makes a single div */
1523 136147098 : long v = 0;
1524 136147098 : if (!r)
1525 : {
1526 22291752 : do { v++; N = q; q = N / p; r = N % p; } while (!r);
1527 14241310 : *n = N;
1528 : }
1529 136147098 : *stop = q <= p; return v;
1530 : }
1531 : /* Assume n > 0. Return v_p(n), set *n := n/p^v_p(n). Set 'stop' if now
1532 : * n < p^2 [implies n prime if no prime < p divides n] */
1533 : long
1534 117961657 : Z_lvalrem_stop(GEN *n, ulong p, int *stop)
1535 : {
1536 : pari_sp av;
1537 : long v;
1538 : ulong r;
1539 : GEN N, q;
1540 :
1541 117961657 : if (lgefint(*n) == 3)
1542 : {
1543 102865413 : r = (*n)[2];
1544 102865413 : v = u_lvalrem_stop(&r, p, stop);
1545 102865413 : if (v) *n = utoipos(r);
1546 102865413 : return v;
1547 : }
1548 15096244 : av = avma; v = 0; q = absdiviu_rem(*n, p, &r);
1549 15096244 : if (r) set_avma(av);
1550 : else
1551 : {
1552 : do {
1553 265866 : v++; N = q;
1554 265866 : if (v == VAL_DC_THRESHOLD)
1555 : {
1556 629 : v += Z_pvalrem_DC(N,sqru(p),&N) << 1;
1557 629 : q = absdiviu_rem(N, p, &r); if (!r) { v++; N = q; }
1558 629 : break;
1559 : }
1560 265237 : q = absdiviu_rem(N, p, &r);
1561 265237 : } while (!r);
1562 224451 : *n = N;
1563 : }
1564 15096244 : *stop = isless_iu(q,p); return v;
1565 : }
1566 :
1567 : /* x is a nonzero integer, |p| > 1 */
1568 : long
1569 69974423 : Z_pvalrem(GEN x, GEN p, GEN *py)
1570 : {
1571 : long vx;
1572 : pari_sp av;
1573 :
1574 69974423 : if (lgefint(p) == 3) return Z_lvalrem(x, uel(p,2), py);
1575 14540555 : if (lgefint(x) == 3) { *py = icopy(x); return 0; }
1576 1575360 : av = avma; vx = 0; (void)new_chunk(lgefint(x));
1577 : for(;;)
1578 21848 : {
1579 1597208 : GEN r, q = dvmdii(x,p,&r);
1580 1597208 : if (r != gen_0) { set_avma(av); *py = icopy(x); return vx; }
1581 21848 : vx++; x = q;
1582 : }
1583 : }
1584 : long
1585 2634708 : u_pvalrem(ulong x, GEN p, ulong *py)
1586 : {
1587 2634708 : if (lgefint(p) == 3) return u_lvalrem(x, uel(p,2), py);
1588 557 : *py = x; return 0;
1589 : }
1590 : long
1591 137019 : u_pval(ulong x, GEN p)
1592 : {
1593 137019 : if (lgefint(p) == 3) return u_lval(x, uel(p,2));
1594 0 : return 0;
1595 : }
1596 : long
1597 50414050 : Z_pval(GEN x, GEN p) {
1598 : long vx;
1599 : pari_sp av;
1600 :
1601 50414050 : if (lgefint(p) == 3) return Z_lval(x, uel(p,2));
1602 32879 : if (lgefint(x) == 3) return 0;
1603 8104 : av = avma; vx = 0;
1604 : for(;;)
1605 25820 : {
1606 33924 : GEN r, q = dvmdii(x,p,&r);
1607 33924 : if (r != gen_0) return gc_long(av,vx);
1608 25820 : vx++; x = q;
1609 : }
1610 : }
1611 :
1612 : /* return v_p(n!) = [n/p] + [n/p^2] + ... */
1613 : long
1614 1995413 : factorial_lval(ulong n, ulong p)
1615 : {
1616 : ulong q, v;
1617 1995413 : if (p == 2) return n - hammingu(n);
1618 1320668 : q = p; v = 0;
1619 1449331 : do { v += n/q; q *= p; } while (n >= q);
1620 1320668 : return (long)v;
1621 : }
1622 :
1623 : /********** Same for "containers" ZX / ZV / ZC **********/
1624 :
1625 : /* If the t_INT q divides the ZX/ZV x, return the quotient. Otherwise NULL.
1626 : * Stack clean; assumes lg(x) > 1 */
1627 : static GEN
1628 6855 : gen_Z_divides(GEN x, GEN q, long imin)
1629 : {
1630 : long i, l;
1631 6855 : GEN y = cgetg_copy(x, &l);
1632 :
1633 6855 : y[1] = x[1]; /* Needed for ZX; no-op if ZV, overwritten in first iteration */
1634 89952 : for (i = imin; i < l; i++)
1635 : {
1636 87235 : GEN r, xi = gel(x,i);
1637 87235 : if (!signe(xi)) { gel(y,i) = xi; continue; }
1638 56650 : gel(y,i) = dvmdii(xi, q, &r);
1639 56650 : if (r != gen_0) { set_avma((pari_sp)(y+l)); return NULL; }
1640 : }
1641 2717 : return y;
1642 : }
1643 : /* If q divides the ZX/ZV x, return the quotient. Otherwise NULL.
1644 : * Stack clean; assumes lg(x) > 1 */
1645 : static GEN
1646 4935 : gen_z_divides(GEN x, ulong q, long imin)
1647 : {
1648 : long i, l;
1649 4935 : GEN y = cgetg_copy(x, &l);
1650 :
1651 4935 : y[1] = x[1]; /* Needed for ZX; no-op if ZV, overwritten in first iteration */
1652 42915 : for (i = imin; i < l; i++)
1653 : {
1654 : ulong r;
1655 41593 : GEN xi = gel(x,i);
1656 41593 : if (!signe(xi)) { gel(y,i) = xi; continue; }
1657 28208 : gel(y,i) = absdiviu_rem(xi, q, &r);
1658 28208 : if (r) { set_avma((pari_sp)(y+l)); return NULL; }
1659 24595 : affectsign_safe(xi, &gel(y,i));
1660 : }
1661 1322 : return y;
1662 : }
1663 :
1664 : /* return v_q(x) and set *py = x / q^v_q(x), using divide & conquer */
1665 : static long
1666 11752 : gen_pvalrem_DC(GEN x, GEN q, GEN *py, long imin)
1667 : {
1668 :
1669 11752 : pari_sp av = avma;
1670 11752 : long v, i, l, lz = LONG_MAX;
1671 11752 : GEN y = cgetg_copy(x, &l);
1672 :
1673 11752 : y[1] = x[1];
1674 136351 : for (i = imin; i < l; i++)
1675 : {
1676 129496 : GEN r, xi = gel(x,i);
1677 129496 : if (!signe(xi)) { gel(y,i) = xi; continue; }
1678 87592 : gel(y,i) = dvmdii(xi, q, &r);
1679 87592 : if (r != gen_0) { *py = x; return gc_long(av,0); }
1680 82695 : lz = minss(lz, lgefint(gel(y,i)));
1681 : }
1682 6855 : if (2 * lgefint(q) <= lz+3) /* avoid squaring if pointless */
1683 6803 : v = gen_pvalrem_DC(y, sqri(q), py, imin) << 1;
1684 : else
1685 52 : { v = 0; *py = y; }
1686 :
1687 6855 : y = gen_Z_divides(*py, q, imin);
1688 6855 : if (!y) return v+1;
1689 2717 : *py = y; return v+2;
1690 : }
1691 :
1692 : static long
1693 816027 : gen_2val(GEN x, long imin)
1694 : {
1695 816027 : long i, lx = lg(x), v = LONG_MAX;
1696 3099359 : for (i = imin; i < lx; i++)
1697 : {
1698 2639347 : GEN c = gel(x,i);
1699 : long w;
1700 2639347 : if (!signe(c)) continue;
1701 2398557 : w = vali(c);
1702 2398557 : if (w < v) { v = w; if (!v) break; }
1703 : }
1704 816027 : return v;
1705 : }
1706 : static long
1707 913697 : gen_lval(GEN x, ulong p, long imin)
1708 : {
1709 : long i, lx, v;
1710 : pari_sp av;
1711 : GEN y;
1712 913697 : if (p == 2) return gen_2val(x, imin);
1713 97670 : av = avma;
1714 97670 : lx = lg(x); y = leafcopy(x);
1715 293227 : for(v = 0;; v++)
1716 1437861 : for (i = imin; i < lx; i++)
1717 : {
1718 : ulong r;
1719 1242304 : gel(y,i) = absdiviu_rem(gel(y,i), p, &r);
1720 1242304 : if (r) return gc_long(av,v);
1721 : }
1722 : }
1723 : long
1724 788796 : ZX_lval(GEN x, ulong p) { return gen_lval(x, p, 2); }
1725 : long
1726 0 : ZV_lval(GEN x, ulong p) { return gen_lval(x, p, 1); }
1727 :
1728 : long
1729 28931 : zx_lval(GEN f, long p)
1730 : {
1731 28931 : long i, l = lg(f), x = LONG_MAX;
1732 30219 : for(i=2; i<l; i++)
1733 : {
1734 : long y;
1735 29477 : if (f[i] == 0) continue;
1736 29428 : y = z_lval(f[i], p);
1737 29428 : if (y < x) { x = y; if (x == 0) return x; }
1738 : }
1739 742 : return x;
1740 : }
1741 :
1742 : static long
1743 135226 : gen_pval(GEN x, GEN p, long imin)
1744 : {
1745 : long i, lx, v;
1746 : pari_sp av;
1747 : GEN y;
1748 135226 : if (lgefint(p) == 3) return gen_lval(x, p[2], imin);
1749 10325 : av = avma;
1750 10325 : lx = lg(x); y = leafcopy(x);
1751 10332 : for(v = 0;; v++)
1752 : {
1753 10332 : if (v == VAL_DC_THRESHOLD)
1754 : {
1755 0 : if (is_pm1(p)) pari_err_DOMAIN("gen_pval", "p", "=", p, p);
1756 0 : v += gen_pvalrem_DC(y, p, &y, imin);
1757 0 : return gc_long(av,v);
1758 : }
1759 :
1760 10381 : for (i = imin; i < lx; i++)
1761 : {
1762 10374 : GEN r; gel(y,i) = dvmdii(gel(y,i), p, &r);
1763 10374 : if (r != gen_0) return gc_long(av,v);
1764 : }
1765 : }
1766 : }
1767 : long
1768 104292 : ZX_pval(GEN x, GEN p) { return gen_pval(x, p, 2); }
1769 : long
1770 30934 : ZV_pval(GEN x, GEN p) { return gen_pval(x, p, 1); }
1771 : /* v = 0 (mod p) */
1772 : int
1773 1414 : ZV_Z_dvd(GEN v, GEN p)
1774 : {
1775 1414 : pari_sp av = avma;
1776 1414 : long i, l = lg(v);
1777 4753 : for (i=1; i<l; i++)
1778 3563 : if (!dvdii(gel(v,i), p)) return gc_int(av, 0);
1779 1190 : return gc_int(av, 1);
1780 : }
1781 :
1782 : static long
1783 4924002 : gen_2valrem(GEN x, GEN *px, long imin)
1784 : {
1785 4924002 : long i, lx = lg(x), v = LONG_MAX;
1786 : GEN z;
1787 14117281 : for (i = imin; i < lx; i++)
1788 : {
1789 12794497 : GEN c = gel(x,i);
1790 : long w;
1791 12794497 : if (!signe(c)) continue;
1792 11897303 : w = vali(c);
1793 11897303 : if (w < v) {
1794 7052142 : v = w;
1795 7052142 : if (!v) { *px = x; return 0; } /* early abort */
1796 : }
1797 : }
1798 1322784 : z = cgetg_copy(x, &lx); z[1] = x[1];
1799 8615669 : for (i=imin; i<lx; i++) gel(z,i) = shifti(gel(x,i), -v);
1800 1322784 : *px = z; return v;
1801 : }
1802 : static long
1803 8489814 : gen_lvalrem(GEN x, ulong p, GEN *px, long imin)
1804 : {
1805 : long i, lx, v;
1806 : GEN y;
1807 8489814 : if (p == 2) return gen_2valrem(x, px, imin);
1808 3565812 : y = cgetg_copy(x, &lx);
1809 3565812 : y[1] = x[1];
1810 3565812 : x = leafcopy(x);
1811 3565812 : for(v = 0;; v++)
1812 : {
1813 4982373 : if (v == VAL_DC_THRESHOLD)
1814 : {
1815 4935 : if (p == 1) pari_err_DOMAIN("gen_lvalrem", "p", "=", gen_1, gen_1);
1816 4935 : v += gen_pvalrem_DC(x, sqru(p), px, imin) << 1;
1817 4935 : x = gen_z_divides(*px, p, imin);
1818 4935 : if (x) { *px = x; v++; }
1819 4935 : return v;
1820 : }
1821 :
1822 15707807 : for (i = imin; i < lx; i++)
1823 : {
1824 14291246 : ulong r; gel(y,i) = absdiviu_rem(gel(x,i), p, &r);
1825 14291246 : if (r) { *px = x; return v; }
1826 10730369 : affectsign_safe(gel(x,i), &gel(y,i));
1827 : }
1828 1416561 : swap(x, y);
1829 : }
1830 : }
1831 : long
1832 721 : ZX_lvalrem(GEN x, ulong p, GEN *px) { return gen_lvalrem(x,p,px, 2); }
1833 : long
1834 0 : ZV_lvalrem(GEN x, ulong p, GEN *px) { return gen_lvalrem(x,p,px, 1); }
1835 :
1836 : static long
1837 8502059 : gen_pvalrem(GEN x, GEN p, GEN *px, long imin)
1838 : {
1839 : long i, lx, v;
1840 : GEN y;
1841 8502059 : if (lgefint(p) == 3) return gen_lvalrem(x, p[2], px, imin);
1842 12966 : y = cgetg_copy(x, &lx);
1843 12966 : y[1] = x[1];
1844 12966 : x = leafcopy(x);
1845 12966 : for(v = 0;; v++)
1846 : {
1847 13755 : if (v == VAL_DC_THRESHOLD)
1848 : {
1849 14 : if (is_pm1(p)) pari_err_DOMAIN("gen_pvalrem", "p", "=", p, p);
1850 14 : return v + gen_pvalrem_DC(x, p, px, imin);
1851 : }
1852 :
1853 22489 : for (i = imin; i < lx; i++)
1854 : {
1855 21700 : GEN r; gel(y,i) = dvmdii(gel(x,i), p, &r);
1856 21700 : if (r != gen_0) { *px = x; return v; }
1857 : }
1858 789 : swap(x, y);
1859 : }
1860 : }
1861 : long
1862 4318531 : ZX_pvalrem(GEN x, GEN p, GEN *px) { return gen_pvalrem(x,p,px, 2); }
1863 : long
1864 4183528 : ZV_pvalrem(GEN x, GEN p, GEN *px) { return gen_pvalrem(x,p,px, 1); }
1865 :
1866 : static long
1867 3801 : ZX_gen_pvalrem(GEN x, GEN p, GEN *px, long imin)
1868 : {
1869 : long i, lx, v;
1870 : GEN y;
1871 3801 : y = cgetg_copy(x, &lx);
1872 3801 : y[1] = x[1];
1873 3801 : x = leafcopy(x);
1874 28126 : for (i = imin; i < lx; i++)
1875 24325 : if (typ(gel(x, i)) != t_INT)
1876 : {
1877 15225 : gel(x, i) = leafcopy(gel(x,i));
1878 15225 : gel(y, i) = leafcopy(gel(x,i));
1879 : }
1880 3801 : for(v = 0;; v++)
1881 : {
1882 : #if 0
1883 : if (v == VAL_DC_THRESHOLD) /* TODO */
1884 : {
1885 : if (is_pm1(p)) pari_err_DOMAIN("ZX_gen_pvalrem", "p", "=", p, p);
1886 : return v + ZX_gen_pvalrem_DC(x, p, px, imin);
1887 : }
1888 : #endif
1889 :
1890 12747 : for (i = imin; i < lx; i++)
1891 : {
1892 11263 : GEN r, xi = gel(x,i);
1893 11263 : if (typ(xi) == t_INT)
1894 : {
1895 3066 : gel(y,i) = dvmdii(xi, p, &r);
1896 5537 : if (r != gen_0) { *px = x; return v; }
1897 : } else
1898 : {
1899 8197 : long j, lxi = lg(xi);
1900 20244 : for(j = 2; j < lxi; j++)
1901 : {
1902 14518 : gmael(y,i,j) = dvmdii(gel(xi,j), p, &r);
1903 14518 : if (r != gen_0) { *px = x; return v; }
1904 : }
1905 : }
1906 : }
1907 1484 : swap(x, y);
1908 : }
1909 : }
1910 :
1911 : long
1912 952 : ZXX_pval(GEN x, GEN p)
1913 : {
1914 952 : pari_sp av = avma;
1915 : GEN px;
1916 952 : return gc_long(av, ZX_gen_pvalrem(x,p,&px, 2));
1917 : }
1918 :
1919 : long
1920 0 : ZXV_pval(GEN x, GEN p)
1921 : {
1922 0 : pari_sp av = avma;
1923 : GEN px;
1924 0 : return gc_long(av, ZX_gen_pvalrem(x,p,&px, 1));
1925 : }
1926 :
1927 : long
1928 2849 : ZXX_pvalrem(GEN x, GEN p, GEN *px) { return ZX_gen_pvalrem(x,p,px, 2); }
1929 : long
1930 0 : ZXV_pvalrem(GEN x, GEN p, GEN *px) { return ZX_gen_pvalrem(x,p,px, 1); }
1931 :
1932 : /*******************************************************************/
1933 : /* */
1934 : /* NEGATION: Create -x */
1935 : /* */
1936 : /*******************************************************************/
1937 :
1938 : GEN
1939 473029126 : gneg(GEN x)
1940 : {
1941 : GEN y;
1942 473029126 : switch(typ(x))
1943 : {
1944 137286518 : case t_INT:
1945 137286518 : return signe(x)? negi(x): gen_0;
1946 247733101 : case t_REAL:
1947 247733101 : return mpneg(x);
1948 :
1949 102754 : case t_INTMOD: y=cgetg(3,t_INTMOD);
1950 102754 : gel(y,1) = icopy(gel(x,1));
1951 102754 : gel(y,2) = signe(gel(x,2))? subii(gel(y,1),gel(x,2)): gen_0;
1952 102754 : break;
1953 :
1954 2617505 : case t_FRAC:
1955 2617505 : y = cgetg(3, t_FRAC);
1956 2617505 : gel(y,1) = negi(gel(x,1));
1957 2617505 : gel(y,2) = icopy(gel(x,2)); break;
1958 :
1959 75286482 : case t_COMPLEX:
1960 75286482 : y=cgetg(3, t_COMPLEX);
1961 75286482 : gel(y,1) = gneg(gel(x,1));
1962 75286482 : gel(y,2) = gneg(gel(x,2));
1963 75286482 : break;
1964 :
1965 250862 : case t_POLMOD:
1966 250862 : retmkpolmod(gneg(gel(x,2)), RgX_copy(gel(x,1)));
1967 :
1968 153370 : case t_RFRAC:
1969 153370 : y = cgetg(3, t_RFRAC);
1970 153370 : gel(y,1) = gneg(gel(x,1));
1971 153370 : gel(y,2) = RgX_copy(gel(x,2)); break;
1972 :
1973 648301 : case t_PADIC:
1974 : {
1975 648301 : GEN u = padic_u(x), pd = padic_pd(x), p = padic_p(x);
1976 648301 : if (!signe(u)) return gcopy(x);
1977 644668 : retmkpadic(subii(pd, u), icopy(p), icopy(pd), valp(x), precp(x));
1978 : }
1979 119 : case t_QUAD:
1980 119 : y=cgetg(4,t_QUAD);
1981 119 : gel(y,1) = ZX_copy(gel(x,1));
1982 119 : gel(y,2) = gneg(gel(x,2));
1983 119 : gel(y,3) = gneg(gel(x,3)); break;
1984 :
1985 2299160 : case t_FFELT: return FF_neg(x);
1986 6214485 : case t_POL: return RgX_neg(x);
1987 16660 : case t_SER: pari_APPLY_ser_normalized(gneg(gel(x,i)));
1988 1813 : case t_VEC: return RgV_neg(x);
1989 431128 : case t_COL: return RgC_neg(x);
1990 1351 : case t_MAT: return RgM_neg(x);
1991 798 : case t_INFINITY: return inf_get_sign(x) == 1? mkmoo(): mkoo();
1992 0 : default:
1993 0 : pari_err_TYPE("gneg",x);
1994 : return NULL; /* LCOV_EXCL_LINE */
1995 : }
1996 78160230 : return y;
1997 : }
1998 :
1999 : GEN
2000 138748057 : gneg_i(GEN x)
2001 : {
2002 : GEN y;
2003 138748057 : switch(typ(x))
2004 : {
2005 64389065 : case t_INT:
2006 64389065 : return signe(x)? negi(x): gen_0;
2007 37512773 : case t_REAL:
2008 37512773 : return mpneg(x);
2009 :
2010 634781 : case t_INTMOD: y=cgetg(3,t_INTMOD);
2011 634781 : gel(y,1) = gel(x,1);
2012 634781 : gel(y,2) = signe(gel(x,2))? subii(gel(y,1),gel(x,2)): gen_0;
2013 634781 : break;
2014 :
2015 5756630 : case t_FRAC:
2016 5756630 : y = cgetg(3, t_FRAC);
2017 5756630 : gel(y,1) = negi(gel(x,1));
2018 5756630 : gel(y,2) = gel(x,2); break;
2019 :
2020 12324287 : case t_COMPLEX:
2021 12324287 : y = cgetg(3, t_COMPLEX);
2022 12324287 : gel(y,1) = gneg_i(gel(x,1));
2023 12324287 : gel(y,2) = gneg_i(gel(x,2)); break;
2024 :
2025 2015218 : case t_PADIC:
2026 : {
2027 2015218 : GEN u = padic_u(x), pd = padic_pd(x), p = padic_p(x);
2028 2015218 : if (!signe(u)) return zeropadic_shallow(p, valp(x));
2029 2013580 : retmkpadic(subii(pd, u), p, pd, valp(x), precp(x));
2030 : }
2031 143436 : case t_POLMOD:
2032 143436 : retmkpolmod(gneg_i(gel(x,2)), RgX_copy(gel(x,1)));
2033 :
2034 85064 : case t_FFELT: return FF_neg_i(x);
2035 :
2036 693 : case t_QUAD: y=cgetg(4,t_QUAD);
2037 693 : gel(y,1) = gel(x,1);
2038 693 : gel(y,2) = gneg_i(gel(x,2));
2039 693 : gel(y,3) = gneg_i(gel(x,3)); break;
2040 :
2041 2730 : case t_VEC:
2042 : case t_COL:
2043 14462 : case t_MAT: pari_APPLY_same(gneg_i(gel(x,i)));
2044 37875357 : case t_POL: pari_APPLY_pol_normalized(gneg_i(gel(x,i)));
2045 4409173 : case t_SER: pari_APPLY_ser_normalized(gneg_i(gel(x,i)));
2046 :
2047 5005055 : case t_RFRAC:
2048 5005055 : y = cgetg(3, t_RFRAC);
2049 5005055 : gel(y,1) = gneg_i(gel(x,1));
2050 5005055 : gel(y,2) = gel(x,2); break;
2051 :
2052 0 : default:
2053 0 : pari_err_TYPE("gneg_i",x);
2054 : return NULL; /* LCOV_EXCL_LINE */
2055 : }
2056 23721446 : return y;
2057 : }
2058 :
2059 : /******************************************************************/
2060 : /* */
2061 : /* ABSOLUTE VALUE */
2062 : /* Create abs(x) if x is integer, real, fraction or complex. */
2063 : /* Error otherwise. */
2064 : /* */
2065 : /******************************************************************/
2066 : static int
2067 0 : is_negative(GEN x) {
2068 0 : switch(typ(x))
2069 : {
2070 0 : case t_INT: case t_REAL:
2071 0 : return (signe(x) < 0);
2072 0 : case t_FRAC:
2073 0 : return (signe(gel(x,1)) < 0);
2074 : }
2075 0 : return 0;
2076 : }
2077 :
2078 : GEN
2079 53743856 : gabs(GEN x, long prec)
2080 : {
2081 : long lx;
2082 : pari_sp av;
2083 : GEN y, N;
2084 :
2085 53743856 : switch(typ(x))
2086 : {
2087 35584068 : case t_INT: case t_REAL:
2088 35584068 : return mpabs(x);
2089 :
2090 13028 : case t_FRAC:
2091 13028 : return absfrac(x);
2092 :
2093 18044169 : case t_COMPLEX:
2094 18044169 : av=avma; N=cxnorm(x);
2095 18044169 : switch(typ(N))
2096 : {
2097 266 : case t_INT:
2098 266 : if (!Z_issquareall(N, &y)) break;
2099 105 : return gc_upto(av, y);
2100 22281 : case t_FRAC: {
2101 : GEN a,b;
2102 36946 : if (!Z_issquareall(gel(N,1), &a)) break;
2103 14665 : if (!Z_issquareall(gel(N,2), &b)) break;
2104 0 : return gc_upto(av, gdiv(a,b));
2105 : }
2106 : }
2107 18044064 : return gc_upto(av, gsqrt(N,prec));
2108 :
2109 21 : case t_QUAD:
2110 21 : av = avma;
2111 21 : return gc_leaf(av, gabs(quadtofp(x, prec), prec));
2112 :
2113 0 : case t_POL:
2114 0 : lx = lg(x); if (lx<=2) return RgX_copy(x);
2115 0 : return is_negative(gel(x,lx-1))? RgX_neg(x): RgX_copy(x);
2116 :
2117 7 : case t_SER:
2118 7 : if (!signe(x)) pari_err_DOMAIN("abs", "argument", "=", gen_0, x);
2119 7 : if (valser(x)) pari_err_DOMAIN("abs", "series valuation", "!=", gen_0, x);
2120 0 : return is_negative(gel(x,2))? gneg(x): gcopy(x);
2121 :
2122 102549 : case t_VEC: case t_COL: case t_MAT:
2123 593784 : pari_APPLY_same(gabs(gel(x,i),prec));
2124 :
2125 14 : case t_INFINITY:
2126 14 : return mkoo();
2127 : }
2128 0 : pari_err_TYPE("gabs",x);
2129 : return NULL; /* LCOV_EXCL_LINE */
2130 : }
2131 :
2132 : GEN
2133 79310 : gmax(GEN x, GEN y) { return gcopy(gmax_shallow(x,y)); }
2134 : GEN
2135 0 : gmaxgs(GEN x, long s) { return (gcmpsg(s,x)>=0)? stoi(s): gcopy(x); }
2136 :
2137 : GEN
2138 12180 : gmin(GEN x, GEN y) { return gcopy(gmin_shallow(x,y)); }
2139 : GEN
2140 0 : gmings(GEN x, long s) { return (gcmpsg(s,x)>0)? gcopy(x): stoi(s); }
2141 :
2142 : long
2143 503449 : vecindexmax(GEN x)
2144 : {
2145 503449 : long lx = lg(x), i0, i;
2146 : GEN s;
2147 :
2148 503449 : if (lx==1) pari_err_DOMAIN("vecindexmax", "empty argument", "=", x,x);
2149 503449 : switch(typ(x))
2150 : {
2151 503449 : case t_VEC: case t_COL:
2152 503449 : s = gel(x,i0=1);
2153 1502952 : for (i=2; i<lx; i++)
2154 999503 : if (gcmp(gel(x,i),s) > 0) s = gel(x,i0=i);
2155 503449 : return i0;
2156 0 : case t_VECSMALL:
2157 0 : return vecsmall_indexmax(x);
2158 0 : default: pari_err_TYPE("vecindexmax",x);
2159 : }
2160 : /* LCOV_EXCL_LINE */
2161 0 : return 0;
2162 : }
2163 : long
2164 183130 : vecindexmin(GEN x)
2165 : {
2166 183130 : long lx = lg(x), i0, i;
2167 : GEN s;
2168 :
2169 183130 : if (lx==1) pari_err_DOMAIN("vecindexmin", "empty argument", "=", x,x);
2170 183130 : switch(typ(x))
2171 : {
2172 183130 : case t_VEC: case t_COL:
2173 183130 : s = gel(x,i0=1);
2174 952909 : for (i=2; i<lx; i++)
2175 769779 : if (gcmp(gel(x,i),s) < 0) s = gel(x,i0=i);
2176 183130 : return i0;
2177 0 : case t_VECSMALL:
2178 0 : return vecsmall_indexmin(x);
2179 0 : default: pari_err_TYPE("vecindexmin",x);
2180 : }
2181 : /* LCOV_EXCL_LINE */
2182 0 : return 0;
2183 : }
2184 :
2185 : GEN
2186 227037 : vecmax0(GEN x, GEN *pi)
2187 : {
2188 227037 : long i, lx = lg(x), tx = typ(x);
2189 227037 : if (!is_matvec_t(tx) && tx != t_VECSMALL
2190 49 : && (tx != t_LIST || list_typ(x) != t_LIST_RAW)) return gcopy(x);
2191 227016 : if (tx == t_LIST)
2192 28 : { if (list_data(x)) { x = list_data(x); lx = lg(x); } else lx = 1; }
2193 227016 : if (lx==1) pari_err_DOMAIN("vecmax", "empty argument", "=", x,x);
2194 226974 : switch(typ(x))
2195 : {
2196 226497 : case t_VEC: case t_COL:
2197 226497 : i = vecindexmax(x); if (pi) *pi = utoipos(i);
2198 226497 : return gcopy(gel(x,i));
2199 456 : case t_MAT: {
2200 456 : long j, i0 = 1, j0 = 1, lx2 = lgcols(x);
2201 : GEN s;
2202 456 : if (lx2 == 1) pari_err_DOMAIN("vecmax", "empty argument", "=", x,x);
2203 449 : s = gcoeff(x,i0,j0); i = 2;
2204 1205 : for (j=1; j<lx; j++,i=1)
2205 : {
2206 756 : GEN c = gel(x,j);
2207 1817 : for (; i<lx2; i++)
2208 1061 : if (gcmp(gel(c,i),s) > 0) { s = gel(c,i); j0=j; i0=i; }
2209 : }
2210 449 : if (pi) *pi = mkvec2(utoipos(i0), utoipos(j0));
2211 449 : return gcopy(s);
2212 : }
2213 21 : case t_VECSMALL:
2214 21 : i = vecsmall_indexmax(x); if (pi) *pi = utoipos(i);
2215 21 : return stoi(x[i]);
2216 : }
2217 : return NULL;/*LCOV_EXCL_LINE*/
2218 : }
2219 : GEN
2220 147521 : vecmin0(GEN x, GEN *pi)
2221 : {
2222 147521 : long i, lx = lg(x), tx = typ(x);
2223 147521 : if (!is_matvec_t(tx) && tx != t_VECSMALL
2224 49 : && (tx != t_LIST || list_typ(x) != t_LIST_RAW)) return gcopy(x);
2225 147500 : if (tx == t_LIST)
2226 28 : { if (list_data(x)) { x = list_data(x); lx = lg(x); } else lx = 1; }
2227 147500 : if (lx==1) pari_err_DOMAIN("vecmin", "empty argument", "=", x,x);
2228 147465 : switch(typ(x))
2229 : {
2230 147423 : case t_VEC: case t_COL:
2231 147423 : i = vecindexmin(x); if (pi) *pi = utoipos(i);
2232 147423 : return gcopy(gel(x,i));
2233 21 : case t_MAT: {
2234 21 : long j, i0 = 1, j0 = 1, lx2 = lgcols(x);
2235 : GEN s;
2236 21 : if (lx2 == 1) pari_err_DOMAIN("vecmin", "empty argument", "=", x,x);
2237 21 : s = gcoeff(x,i0,j0); i = 2;
2238 63 : for (j=1; j<lx; j++,i=1)
2239 : {
2240 42 : GEN c = gel(x,j);
2241 105 : for (; i<lx2; i++)
2242 63 : if (gcmp(gel(c,i),s) < 0) { s = gel(c,i); j0=j; i0=i; }
2243 : }
2244 21 : if (pi) *pi = mkvec2(utoipos(i0), utoipos(j0));
2245 21 : return gcopy(s);
2246 : }
2247 21 : case t_VECSMALL:
2248 21 : i = vecsmall_indexmin(x); if (pi) *pi = utoipos(i);
2249 21 : return stoi(x[i]);
2250 : }
2251 : return NULL;/*LCOV_EXCL_LINE*/
2252 : }
2253 :
2254 : GEN
2255 66207 : vecmax(GEN x) { return vecmax0(x, NULL); }
2256 : GEN
2257 35353 : vecmin(GEN x) { return vecmin0(x, NULL); }
2258 :
2259 : /*******************************************************************/
2260 : /* */
2261 : /* GENERIC AFFECTATION */
2262 : /* Affect the content of x to y, whenever possible */
2263 : /* */
2264 : /*******************************************************************/
2265 : /* x PADIC, Y INT, return lift(x * Mod(1,Y)) */
2266 : GEN
2267 4550 : padic_to_Fp(GEN x, GEN Y) {
2268 4550 : pari_sp av = avma;
2269 4550 : GEN p = padic_p(x), z;
2270 4550 : long vy, vx = valp(x);
2271 4550 : if (!signe(Y)) pari_err_INV("padic_to_Fp",Y);
2272 4550 : vy = Z_pvalrem(Y,p, &z);
2273 4550 : if (vx < 0 || !gequal1(z)) pari_err_OP("",x, mkintmod(gen_1,Y));
2274 4529 : if (vx >= vy) { set_avma(av); return gen_0; }
2275 4172 : z = padic_u(x);
2276 4172 : if (!signe(z) || vy > vx + precp(x)) pari_err_OP("",x, mkintmod(gen_1,Y));
2277 4172 : if (vx) z = mulii(z, powiu(p,vx));
2278 4172 : return gc_INT(av, remii(z, Y));
2279 : }
2280 : ulong
2281 421719 : padic_to_Fl(GEN x, ulong Y) {
2282 421719 : GEN p = padic_p(x);
2283 : ulong u, z;
2284 421719 : long vy, vx = valp(x);
2285 421719 : vy = u_pvalrem(Y,p, &u);
2286 421719 : if (vx < 0 || u != 1) pari_err_OP("",x, mkintmodu(1,Y));
2287 : /* Y = p^vy */
2288 421719 : if (vx >= vy) return 0;
2289 376009 : z = umodiu(padic_u(x), Y);
2290 376009 : if (!z || vy > vx + precp(x)) pari_err_OP("",x, mkintmodu(1,Y));
2291 376009 : if (vx) {
2292 0 : ulong pp = p[2];
2293 0 : z = Fl_mul(z, upowuu(pp,vx), Y); /* p^vx < p^vy = Y */
2294 : }
2295 376009 : return z;
2296 : }
2297 :
2298 : /* y a t_COMPLEX of t_REAL from cgetc; x a "complex number" (INT, FRAC, REAL,
2299 : * COMPLEX), raise an exception on exotic types (t_QUAD, etc) */
2300 : void
2301 2590 : affgc(GEN x, GEN y)
2302 : {
2303 2590 : switch(typ(x))
2304 : {
2305 0 : case t_INT:
2306 0 : affir(x, gel(y,1));
2307 0 : affir(gen_0,gel(y,2)); break;
2308 :
2309 910 : case t_REAL:
2310 910 : affrr(x,gel(y,1));
2311 910 : affir(gen_0,gel(y,2)); break;
2312 :
2313 1680 : case t_COMPLEX:
2314 1680 : affgr(gel(x,1),gel(y,1));
2315 1680 : affgr(gel(x,2),gel(y,2)); break;
2316 :
2317 0 : case t_FRAC:
2318 0 : rdiviiz(gel(x,1),gel(x,2), gel(y,1));
2319 0 : affir(gen_0,gel(y,2)); break;
2320 :
2321 0 : default: pari_err_TYPE2("=",x,y);
2322 : }
2323 2590 : }
2324 :
2325 : /*******************************************************************/
2326 : /* */
2327 : /* CONVERSION QUAD --> REAL, COMPLEX OR P-ADIC */
2328 : /* */
2329 : /*******************************************************************/
2330 : GEN
2331 252 : quadtofp(GEN x, long prec)
2332 : {
2333 252 : GEN b, D, z, u = gel(x,2), v = gel(x,3);
2334 : pari_sp av;
2335 252 : if (prec < LOWDEFAULTPREC) prec = LOWDEFAULTPREC;
2336 252 : if (isintzero(v)) return cxcompotor(u, prec);
2337 252 : av = avma; D = quad_disc(x); b = gel(gel(x,1),3); /* 0 or -1 */
2338 : /* u + v (-b + sqrt(D)) / 2 */
2339 252 : if (!signe(b)) b = NULL;
2340 252 : if (b) u = gadd(gmul2n(u,1), v);
2341 252 : z = sqrtr_abs(itor(D, prec));
2342 252 : if (!b) shiftr_inplace(z, -1);
2343 252 : z = gmul(v, z);
2344 252 : if (signe(D) < 0)
2345 : {
2346 35 : z = mkcomplex(cxcompotor(u, prec), z);
2347 35 : if (!b) return gc_GEN(av, z);
2348 0 : z = gmul2n(z, -1);
2349 : }
2350 : else
2351 : { /* if (b) x ~ (u + z) / 2 and quadnorm(x) ~ (u^2 - z^2) / 4
2352 : * else x ~ u + z and quadnorm(x) ~ u^2 - z^2 */
2353 217 : long s = gsigne(u);
2354 217 : if (s == -gsigne(v)) /* conjugate expression avoids cancellation */
2355 : {
2356 14 : z = gdiv(quadnorm(x), gsub(u, z));
2357 14 : if (b) shiftr_inplace(z, 1);
2358 : }
2359 : else
2360 : {
2361 203 : if (s) z = gadd(u, z);
2362 203 : if (b) shiftr_inplace(z, -1);
2363 : }
2364 : }
2365 217 : return gc_upto(av, z);
2366 : }
2367 :
2368 : static GEN
2369 42 : qtop(GEN x, GEN p, long d)
2370 : {
2371 42 : GEN z, D, P, b, u = gel(x,2), v = gel(x,3);
2372 : pari_sp av;
2373 42 : if (gequal0(v)) return cvtop(u, p, d);
2374 28 : P = gel(x,1);
2375 28 : b = gel(P,3);
2376 28 : av = avma; D = quad_disc(x);
2377 28 : if (absequaliu(p,2)) d += 2;
2378 28 : z = Qp_sqrt(cvtop(D,p,d));
2379 28 : if (!z) pari_err_SQRTN("Qp_sqrt",D);
2380 14 : z = gmul2n(gsub(z, b), -1);
2381 :
2382 14 : z = gadd(u, gmul(v, z));
2383 14 : if (typ(z) != t_PADIC) /* t_INTMOD for t_QUAD of t_INTMODs... */
2384 0 : z = cvtop(z, p, d);
2385 14 : return gc_upto(av, z);
2386 : }
2387 : static GEN
2388 14 : ctop(GEN x, GEN p, long d)
2389 : {
2390 14 : pari_sp av = avma;
2391 14 : GEN z, u = gel(x,1), v = gel(x,2);
2392 14 : if (isrationalzero(v)) return cvtop(u, p, d);
2393 14 : z = Qp_sqrt(cvtop(gen_m1, p, d - gvaluation(v, p))); /* = I */
2394 14 : if (!z) pari_err_SQRTN("Qp_sqrt",gen_m1);
2395 :
2396 14 : z = gadd(u, gmul(v, z));
2397 14 : if (typ(z) != t_PADIC) /* t_INTMOD for t_COMPLEX of t_INTMODs... */
2398 0 : z = cvtop(z, p, d);
2399 14 : return gc_upto(av, z);
2400 : }
2401 :
2402 : /* cvtop2(stoi(s), y) */
2403 : GEN
2404 399 : cvstop2(long s, GEN y)
2405 : {
2406 399 : GEN p = padic_p(y), pd = padic_pd(y), u = padic_u(y);
2407 399 : long v, d = signe(u)? precp(y): 0;
2408 399 : if (!s) return zeropadic_shallow(p, d);
2409 399 : v = z_pvalrem(s, p, &s);
2410 399 : if (d <= 0) return zeropadic_shallow(p, v);
2411 399 : retmkpadic(modsi(s, pd), p, pd, v, d);
2412 : }
2413 :
2414 : static GEN
2415 17618519 : itop2_coprime(GEN x, GEN y, long v, long d)
2416 : {
2417 17618519 : GEN p = padic_p(y), pd = padic_pd(y);
2418 17618519 : retmkpadic(modii(x, pd), p, pd, v, d);
2419 : }
2420 : /* cvtop(x, gel(y,2), precp(y)), shallow */
2421 : GEN
2422 17625449 : cvtop2(GEN x, GEN y)
2423 : {
2424 17625449 : GEN p = padic_p(y), u;
2425 17625449 : long v, d = signe(padic_u(y))? precp(y): 0;
2426 17625449 : switch(typ(x))
2427 : {
2428 14826689 : case t_INT:
2429 14826689 : if (!signe(x)) return zeropadic_shallow(p, d);
2430 14826689 : if (d <= 0) return zeropadic_shallow(p, Z_pval(x,p));
2431 14822251 : v = Z_pvalrem(x, p, &x); return itop2_coprime(x, y, v, d);
2432 :
2433 0 : case t_INTMOD:
2434 0 : v = Z_pval(gel(x,1),p); if (v > d) v = d;
2435 0 : return cvtop(gel(x,2), p, v);
2436 :
2437 2797514 : case t_FRAC:
2438 : {
2439 : GEN num, den;
2440 2797514 : if (d <= 0) return zeropadic_shallow(p, Q_pval(x,p));
2441 2796233 : num = gel(x,1); v = Z_pvalrem(num, p, &num);
2442 2796233 : den = gel(x,2); if (!v) v = -Z_pvalrem(den, p, &den);
2443 2796233 : if (!is_pm1(den)) num = mulii(num, Fp_inv(den, padic_pd(y)));
2444 2796233 : return itop2_coprime(num, y, v, d);
2445 : }
2446 7 : case t_COMPLEX: return ctop(x, p, d);
2447 42 : case t_QUAD: return qtop(x, p, d);
2448 1197 : case t_PADIC:
2449 1197 : u = padic_u(x);
2450 1197 : if (!signe(u)) return zeropadic_shallow(p, d);
2451 1197 : if (precp(x) <= d) return x;
2452 35 : return itop2_coprime(u, y, valp(x), d); /* reduce accuracy */
2453 : }
2454 0 : pari_err_TYPE("cvtop2",x);
2455 : return NULL; /* LCOV_EXCL_LINE */
2456 : }
2457 :
2458 : static GEN
2459 163005 : _Fp_div(GEN n, GEN d, GEN q)
2460 163005 : { return equali1(d)? modii(n, q): Fp_div(n, d, q); }
2461 :
2462 : /* assume is_const_t(tx) */
2463 : GEN
2464 600636 : cvtop(GEN x, GEN p, long d)
2465 : {
2466 : GEN u;
2467 : long v;
2468 :
2469 600636 : if (typ(p) != t_INT) pari_err_TYPE("cvtop",p);
2470 600636 : switch(typ(x))
2471 : {
2472 268304 : case t_INT:
2473 268304 : if (!signe(x)) return zeropadic(p, d);
2474 267086 : if (d <= 0) return zeropadic(p, Z_pval(x,p));
2475 267030 : v = Z_pvalrem(x, p, &x); /* not memory-clean */
2476 267030 : retmkpadic_i(modii(x, _pd), icopy(p), powiu(p,d), v, d);
2477 :
2478 28 : case t_INTMOD:
2479 28 : v = Z_pval(gel(x,1),p); if (v > d) v = d;
2480 28 : return cvtop(gel(x,2), p, v);
2481 :
2482 163019 : case t_FRAC:
2483 : {
2484 : GEN num, den;
2485 163019 : if (d <= 0) return zeropadic(p, Q_pval(x,p));
2486 163005 : num = gel(x,1); v = Z_pvalrem(num, p, &num); /* not memory-clean */
2487 163005 : den = gel(x,2); if (!v) v = -Z_pvalrem(den, p, &den);
2488 163005 : retmkpadic_i(_Fp_div(num, den, _pd), icopy(p), powiu(p,d), v, d);
2489 : }
2490 7 : case t_COMPLEX: return ctop(x, p, d);
2491 169278 : case t_PADIC:
2492 169278 : p = padic_p(x); /* override */
2493 169278 : u = padic_u(x);
2494 169278 : if (!signe(u)) return zeropadic(p, d);
2495 169194 : retmkpadic_i(modii(u, _pd), icopy(p), powiu(p,d), valp(x), d);
2496 :
2497 0 : case t_QUAD: return qtop(x, p, d);
2498 : }
2499 0 : pari_err_TYPE("cvtop",x);
2500 : return NULL; /* LCOV_EXCL_LINE */
2501 : }
2502 :
2503 : GEN
2504 126 : gcvtop(GEN x, GEN p, long r)
2505 : {
2506 126 : switch(typ(x))
2507 : {
2508 63 : case t_POL: pari_APPLY_pol_normalized(gcvtop(gel(x,i),p,r));
2509 35 : case t_SER: pari_APPLY_ser_normalized(gcvtop(gel(x,i),p,r));
2510 0 : case t_POLMOD: case t_RFRAC: case t_VEC: case t_COL: case t_MAT:
2511 0 : pari_APPLY_same(gcvtop(gel(x,i),p,r));
2512 : }
2513 98 : return cvtop(x,p,r);
2514 : }
2515 :
2516 : long
2517 890037500 : gexpo_safe(GEN x)
2518 : {
2519 890037500 : long tx = typ(x), lx, e, f, i;
2520 :
2521 890037500 : switch(tx)
2522 : {
2523 200001587 : case t_INT:
2524 200001587 : return expi(x);
2525 :
2526 646152 : case t_FRAC:
2527 646152 : return expi(gel(x,1)) - expi(gel(x,2));
2528 :
2529 478490917 : case t_REAL:
2530 478490917 : return expo(x);
2531 :
2532 88548937 : case t_COMPLEX:
2533 88548937 : e = gexpo(gel(x,1));
2534 88548937 : f = gexpo(gel(x,2)); return maxss(e, f);
2535 :
2536 91 : case t_QUAD: {
2537 91 : GEN p = gel(x,1); /* mod = X^2 + {0,1}* X - {D/4, (1-D)/4})*/
2538 91 : long d = 1 + expi(gel(p,2))/2; /* ~ expo(sqrt(D)) */
2539 91 : e = gexpo(gel(x,2));
2540 91 : f = gexpo(gel(x,3)) + d; return maxss(e, f);
2541 : }
2542 84411131 : case t_POL: case t_SER:
2543 84411131 : lx = lg(x); f = -(long)HIGHEXPOBIT;
2544 327534407 : for (i=2; i<lx; i++) { e=gexpo(gel(x,i)); if (e>f) f=e; }
2545 84411131 : return f;
2546 37938629 : case t_VEC: case t_COL: case t_MAT:
2547 37938629 : lx = lg(x); f = -(long)HIGHEXPOBIT;
2548 262116049 : for (i=1; i<lx; i++) { e=gexpo(gel(x,i)); if (e>f) f=e; }
2549 37938629 : return f;
2550 : }
2551 56 : return -1-(long)HIGHEXPOBIT;
2552 : }
2553 : long
2554 889750466 : gexpo(GEN x)
2555 : {
2556 889750466 : long e = gexpo_safe(x);
2557 889750466 : if (e < -(long)HIGHEXPOBIT) pari_err_TYPE("gexpo",x);
2558 889750459 : return e;
2559 : }
2560 : GEN
2561 90285 : gpexponent(GEN x)
2562 : {
2563 90285 : long e = gexpo(x);
2564 90285 : return e == -(long)HIGHEXPOBIT? mkmoo(): stoi(e);
2565 : }
2566 :
2567 : long
2568 7 : sizedigit(GEN x)
2569 : {
2570 7 : return gequal0(x)? 0: (long) ((gexpo(x)+1) * LOG10_2) + 1;
2571 : }
2572 :
2573 : /* normalize series. avma is not updated */
2574 : GEN
2575 15686279 : normalizeser(GEN x)
2576 : {
2577 15686279 : long i, lx = lg(x), vx=varn(x), vp=valser(x);
2578 : GEN y, z;
2579 :
2580 15686279 : if (lx == 2) { setsigne(x,0); return x; }
2581 15685915 : if (lx == 3) {
2582 118989 : z = gel(x,2);
2583 118989 : if (!gequal0(z)) { setsigne(x,1); return x; }
2584 23779 : if (isrationalzero(z)) return zeroser(vx,vp+1);
2585 4256 : if (isexactzero(z)) {
2586 : /* dangerous case: already normalized ? */
2587 91 : if (!signe(x)) return x;
2588 35 : setvalser(x,vp+1); /* no: normalize */
2589 : }
2590 4200 : setsigne(x,0); return x;
2591 : }
2592 15868280 : for (i=2; i<lx; i++)
2593 15821534 : if (! isrationalzero(gel(x,i))) break;
2594 15566926 : if (i == lx) return zeroser(vx,lx-2+vp);
2595 15520180 : z = gel(x,i);
2596 15523995 : while (i<lx && isexactzero(gel(x,i))) i++;
2597 15520180 : if (i == lx)
2598 : {
2599 280 : i -= 3; y = x + i;
2600 280 : stackdummy((pari_sp)y, (pari_sp)x);
2601 280 : gel(y,2) = z;
2602 280 : y[1] = evalsigne(0) | evalvalser(lx-2+vp) | evalvarn(vx);
2603 280 : y[0] = evaltyp(t_SER) | _evallg(3);
2604 280 : return y;
2605 : }
2606 :
2607 15519900 : i -= 2; y = x + i; lx -= i;
2608 15519900 : y[1] = evalsigne(1) | evalvalser(vp+i) | evalvarn(vx);
2609 15519900 : y[0] = evaltyp(t_SER) | _evallg(lx);
2610 :
2611 15519900 : stackdummy((pari_sp)y, (pari_sp)x);
2612 15555546 : for (i = 2; i < lx; i++)
2613 15554279 : if (!gequal0(gel(y, i))) return y;
2614 1267 : setsigne(y, 0); return y;
2615 : }
2616 :
2617 : GEN
2618 0 : normalizepol_approx(GEN x, long lx)
2619 : {
2620 : long i;
2621 0 : for (i = lx-1; i>1; i--)
2622 0 : if (! gequal0(gel(x,i))) break;
2623 0 : stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + i+1));
2624 0 : setlg(x, i+1); setsigne(x, i!=1); return x;
2625 : }
2626 :
2627 : GEN
2628 946183589 : normalizepol_lg(GEN x, long lx)
2629 : {
2630 946183589 : long i, LX = 0;
2631 946183589 : GEN KEEP = NULL;
2632 :
2633 1294541309 : for (i = lx-1; i>1; i--)
2634 : {
2635 1137058764 : GEN z = gel(x,i);
2636 1137058764 : if (! gequal0(z) ) {
2637 788701044 : if (!LX) LX = i+1;
2638 788701044 : stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + LX));
2639 788701044 : x[0] = evaltyp(t_POL) | _evallg(LX);
2640 788701044 : setsigne(x,1); return x;
2641 348357720 : } else if (!isexactzero(z)) {
2642 995968 : if (!LX) LX = i+1; /* to be kept as leading coeff */
2643 347361752 : } else if (!isrationalzero(z))
2644 718617 : KEEP = z; /* to be kept iff all other coeffs are exact 0s */
2645 : }
2646 157482545 : if (!LX) {
2647 157072342 : if (KEEP) { /* e.g. Pol(Mod(0,2)) */
2648 151278 : gel(x,2) = KEEP;
2649 151278 : LX = 3;
2650 : } else
2651 156921064 : LX = 2; /* Pol(0) */
2652 : }
2653 157482545 : stackdummy((pari_sp)(x + lg(x)), (pari_sp)(x + LX));
2654 157482545 : x[0] = evaltyp(t_POL) | _evallg(LX);
2655 157482545 : setsigne(x,0); return x;
2656 : }
2657 :
2658 : /* normalize polynomial x in place */
2659 : GEN
2660 99744606 : normalizepol(GEN x)
2661 : {
2662 99744606 : return normalizepol_lg(x, lg(x));
2663 : }
2664 :
2665 : int
2666 84282121 : gsigne(GEN x)
2667 : {
2668 84282121 : switch(typ(x))
2669 : {
2670 83901329 : case t_INT: case t_REAL: return signe(x);
2671 380092 : case t_FRAC: return signe(gel(x,1));
2672 686 : case t_QUAD:
2673 : {
2674 686 : pari_sp av = avma;
2675 686 : GEN T = gel(x,1), a = gel(x,2), b = gel(x,3);
2676 : long sa, sb;
2677 686 : if (signe(gel(T,2)) > 0) break;
2678 672 : a = gmul2n(a,1);
2679 672 : if (signe(gel(T,3))) a = gadd(a,b);
2680 : /* a + b sqrt(D) > 0 ? */
2681 672 : sa = gsigne(a);
2682 672 : sb = gsigne(b); if (sa == sb) return gc_int(av,sa);
2683 280 : if (sa == 0) return gc_int(av,sb);
2684 266 : if (sb == 0) return gc_int(av,sa);
2685 : /* different signs, take conjugate expression */
2686 259 : sb = gsigne(gsub(gsqr(a), gmul(quad_disc(x), gsqr(b))));
2687 259 : return gc_int(av, sb*sa);
2688 : }
2689 14 : case t_INFINITY: return inf_get_sign(x);
2690 : }
2691 14 : pari_err_TYPE("gsigne",x);
2692 : return 0; /* LCOV_EXCL_LINE */
2693 : }
2694 :
2695 : /*******************************************************************/
2696 : /* */
2697 : /* LISTS */
2698 : /* */
2699 : /*******************************************************************/
2700 : /* make sure L can hold l elements, at least doubling the previous max number
2701 : * of components. */
2702 : static void
2703 810971 : ensure_nb(GEN L, long l)
2704 : {
2705 810971 : long nmax = list_nmax(L), i, lw;
2706 : GEN v, w;
2707 810971 : if (l <= nmax) return;
2708 1246 : if (nmax)
2709 : {
2710 490 : nmax <<= 1;
2711 490 : if (l > nmax) nmax = l;
2712 490 : w = list_data(L); lw = lg(w);
2713 490 : v = newblock(nmax+1);
2714 490 : v[0] = w[0];
2715 1070958 : for (i=1; i < lw; i++) gel(v,i) = gel(w, i);
2716 490 : killblock(w);
2717 : }
2718 : else /* unallocated */
2719 : {
2720 756 : nmax = 32;
2721 756 : if (list_data(L))
2722 0 : pari_err(e_MISC, "store list in variable before appending elements");
2723 756 : v = newblock(nmax+1);
2724 756 : v[0] = evaltyp(t_VEC) | _evallg(1);
2725 : }
2726 1246 : list_data(L) = v;
2727 1246 : L[1] = evaltyp(list_typ(L))|evallg(nmax);
2728 : }
2729 :
2730 : GEN
2731 6926 : mklist_typ(long t)
2732 : {
2733 6926 : GEN L = cgetg(3,t_LIST);
2734 6926 : L[1] = evaltyp(t);
2735 6926 : list_data(L) = NULL; return L;
2736 : }
2737 :
2738 : GEN
2739 6870 : mklist(void)
2740 : {
2741 6870 : return mklist_typ(t_LIST_RAW);
2742 : }
2743 :
2744 : GEN
2745 49 : mkmap(void)
2746 : {
2747 49 : return mklist_typ(t_LIST_MAP);
2748 : }
2749 :
2750 : /* return a list with single element x, allocated on stack */
2751 : GEN
2752 70 : mklistcopy(GEN x)
2753 : {
2754 70 : GEN y = mklist();
2755 70 : list_data(y) = mkveccopy(x);
2756 70 : return y;
2757 : }
2758 :
2759 : GEN
2760 776125 : listput(GEN L, GEN x, long index)
2761 : {
2762 : long l;
2763 : GEN z;
2764 :
2765 776125 : if (index < 0) pari_err_COMPONENT("listput", "<", gen_0, stoi(index));
2766 776118 : z = list_data(L);
2767 776118 : l = z? lg(z): 1;
2768 :
2769 776118 : x = gclone(x);
2770 776118 : if (!index || index >= l)
2771 : {
2772 775964 : ensure_nb(L, l);
2773 775964 : z = list_data(L); /* it may change ! */
2774 775964 : index = l;
2775 775964 : l++;
2776 : } else
2777 154 : gunclone_deep( gel(z, index) );
2778 776118 : gel(z,index) = x;
2779 776118 : z[0] = evaltyp(t_VEC) | evallg(l); /*must be after gel(z,index) is set*/
2780 776118 : return gel(z,index);
2781 : }
2782 :
2783 : GEN
2784 724955 : listput0(GEN L, GEN x, long index)
2785 : {
2786 724955 : if (typ(L) != t_LIST || list_typ(L) != t_LIST_RAW)
2787 14 : pari_err_TYPE("listput",L);
2788 724941 : (void) listput(L, x, index);
2789 724934 : return x;
2790 : }
2791 :
2792 : GEN
2793 35014 : listinsert(GEN L, GEN x, long index)
2794 : {
2795 : long l, i;
2796 : GEN z;
2797 :
2798 35014 : z = list_data(L); l = z? lg(z): 1;
2799 35014 : if (index <= 0) pari_err_COMPONENT("listinsert", "<=", gen_0, stoi(index));
2800 35007 : if (index > l) index = l;
2801 35007 : ensure_nb(L, l);
2802 35007 : BLOCK_SIGINT_START
2803 35007 : z = list_data(L);
2804 87552507 : for (i=l; i > index; i--) gel(z,i) = gel(z,i-1);
2805 35007 : z[0] = evaltyp(t_VEC) | evallg(l+1);
2806 35007 : gel(z,index) = gclone(x);
2807 35007 : BLOCK_SIGINT_END
2808 35007 : return gel(z,index);
2809 : }
2810 :
2811 : GEN
2812 35028 : listinsert0(GEN L, GEN x, long index)
2813 : {
2814 35028 : if (typ(L) != t_LIST || list_typ(L) != t_LIST_RAW)
2815 14 : pari_err_TYPE("listinsert",L);
2816 35014 : (void) listinsert(L, x, index);
2817 35007 : return x;
2818 : }
2819 :
2820 : void
2821 21917 : listpop(GEN L, long index)
2822 : {
2823 : long l, i;
2824 : GEN z;
2825 :
2826 21917 : if (typ(L) != t_LIST) pari_err_TYPE("listinsert",L);
2827 21917 : if (index < 0) pari_err_COMPONENT("listpop", "<", gen_0, stoi(index));
2828 21917 : z = list_data(L);
2829 21917 : if (!z || (l = lg(z)-1) == 0) return;
2830 :
2831 21903 : if (!index || index > l) index = l;
2832 21903 : BLOCK_SIGINT_START
2833 21903 : gunclone_deep( gel(z, index) );
2834 21903 : z[0] = evaltyp(t_VEC) | _evallg(l);
2835 21910 : for (i=index; i < l; i++) z[i] = z[i+1];
2836 21903 : BLOCK_SIGINT_END
2837 : }
2838 :
2839 : void
2840 56 : listpop0(GEN L, long index)
2841 : {
2842 56 : if (typ(L) != t_LIST || list_typ(L) != t_LIST_RAW)
2843 14 : pari_err_TYPE("listpop",L);
2844 42 : listpop(L, index);
2845 42 : }
2846 :
2847 : /* return a copy fully allocated on stack. gclone from changevalue is
2848 : * supposed to malloc() it */
2849 : GEN
2850 5972 : gtolist(GEN x)
2851 : {
2852 : GEN y;
2853 :
2854 5972 : if (!x) return mklist();
2855 391 : switch(typ(x))
2856 : {
2857 314 : case t_VEC: case t_COL:
2858 314 : y = mklist();
2859 314 : if (lg(x) == 1) return y;
2860 293 : list_data(y) = gcopy(x);
2861 293 : settyp(list_data(y), t_VEC);
2862 293 : return y;
2863 7 : case t_LIST:
2864 7 : y = mklist();
2865 7 : list_data(y) = list_data(x)? gcopy(list_data(x)): NULL;
2866 7 : return y;
2867 70 : default:
2868 70 : return mklistcopy(x);
2869 : }
2870 : }
2871 :
2872 : void
2873 21 : listsort(GEN L, long flag)
2874 : {
2875 : long i, l;
2876 21 : pari_sp av = avma;
2877 : GEN perm, v, vnew;
2878 :
2879 21 : if (typ(L) != t_LIST) pari_err_TYPE("listsort",L);
2880 21 : v = list_data(L); l = v? lg(v): 1;
2881 21 : if (l < 3) return;
2882 21 : if (flag)
2883 : {
2884 : long lnew;
2885 14 : perm = gen_indexsort_uniq(L, (void*)&cmp_universal, cmp_nodata);
2886 14 : lnew = lg(perm); /* may have changed since 'uniq' */
2887 14 : vnew = cgetg(lnew,t_VEC);
2888 56 : for (i=1; i<lnew; i++) {
2889 42 : long c = perm[i];
2890 42 : gel(vnew,i) = gel(v,c);
2891 42 : gel(v,c) = NULL;
2892 : }
2893 14 : if (l != lnew) { /* was shortened */
2894 105 : for (i=1; i<l; i++)
2895 91 : if (gel(v,i)) gunclone_deep(gel(v,i));
2896 14 : l = lnew;
2897 : }
2898 : }
2899 : else
2900 : {
2901 7 : perm = gen_indexsort(L, (void*)&cmp_universal, cmp_nodata);
2902 7 : vnew = cgetg(l,t_VEC);
2903 63 : for (i=1; i<l; i++) gel(vnew,i) = gel(v,perm[i]);
2904 : }
2905 119 : for (i=1; i<l; i++) gel(v,i) = gel(vnew,i);
2906 21 : v[0] = vnew[0]; set_avma(av);
2907 : }
|
