Line data Source code
1 : #line 2 "../src/kernel/none/level1.h"
2 : /* Copyright (C) 2000 The PARI group.
3 :
4 : This file is part of the PARI/GP package.
5 :
6 : PARI/GP is free software; you can redistribute it and/or modify it under the
7 : terms of the GNU General Public License as published by the Free Software
8 : Foundation; either version 2 of the License, or (at your option) any later
9 : version. It is distributed in the hope that it will be useful, but WITHOUT
10 : ANY WARRANTY WHATSOEVER.
11 :
12 : Check the License for details. You should have received a copy of it, along
13 : with the package; see the file 'COPYING'. If not, write to the Free Software
14 : Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. */
15 :
16 : /* This file defines "level 1" kernel functions.
17 : * These functions can be inline; they are also defined externally in
18 : * mpinl.c, which includes this file and never needs to be changed */
19 :
20 : INLINE long
21 79500425064 : evallg(long x)
22 : {
23 79500425064 : if (x & ~LGBITS) pari_err_OVERFLOW("lg()");
24 79504069278 : return _evallg(x);
25 : }
26 : INLINE long
27 46093089 : evalvalp(long x)
28 : {
29 46093089 : long v = _evalvalp(x);
30 46093089 : if (v & ~VALPBITS) pari_err_OVERFLOW("valp()");
31 46092650 : return v;
32 : }
33 : INLINE long
34 21139811 : evalvalser(long x)
35 : {
36 21139811 : long v = _evalvalser(x);
37 21139811 : if (v & ~VALSERBITS) pari_err_OVERFLOW("valser()");
38 21139811 : return v;
39 : }
40 : INLINE long
41 10132122557 : evalexpo(long x)
42 : {
43 10132122557 : long v = _evalexpo(x);
44 10132122557 : if (v & ~EXPOBITS) pari_err_OVERFLOW("expo()");
45 10127515014 : return v;
46 : }
47 : INLINE long
48 25050151 : evalprecp(long x)
49 : {
50 25050151 : long v = _evalprecp(x);
51 25050151 : if (x & ~((1UL<<(BITS_IN_LONG-VALPnumBITS))-1)) pari_err_OVERFLOW("precp()");
52 25050214 : return v;
53 : }
54 :
55 : INLINE int
56 162907251 : varncmp(long x, long y)
57 : {
58 162907251 : if (varpriority[x] < varpriority[y]) return 1;
59 125919927 : if (varpriority[x] > varpriority[y]) return -1;
60 68546510 : return 0;
61 : }
62 : INLINE long
63 0 : varnmin(long x, long y)
64 0 : { return (varpriority[x] <= varpriority[y])? x: y; }
65 : INLINE long
66 203 : varnmax(long x, long y)
67 203 : { return (varpriority[x] >= varpriority[y])? x: y; }
68 :
69 : /* Inhibit some area gerepile-wise: declare it to be a non recursive
70 : * type, of length l. Thus gerepile won't inspect the zone, just copy it.
71 : * For the following situation:
72 : * z = cgetg(t,a); av = avma; garbage(); ltop = avma;
73 : * for (i=1; i<HUGE; i++) gel(z,i) = blah();
74 : * stackdummy(av,ltop);
75 : * loses (av-ltop) words but save a costly gerepile. */
76 : INLINE void
77 3081394037 : stackdummy(pari_sp av, pari_sp ltop) {
78 3081394037 : long l = ((GEN)av) - ((GEN)ltop);
79 3081394037 : if (l > 0) {
80 1045145375 : GEN z = (GEN)ltop;
81 1045145375 : z[0] = evaltyp(t_VECSMALL) | evallg(l);
82 : #ifdef DEBUG
83 : { long i; for (i = 1; i < l; i++) z[i] = 0; }
84 : #endif
85 : }
86 3081203311 : }
87 : INLINE void
88 89140888 : fixlg(GEN x, long ly) {
89 89140888 : long lx = lg(x), l = lx - ly;
90 89140888 : if (l > 0)
91 : { /* stackdummy(x+lx, x+ly) */
92 46920486 : GEN z = x + ly;
93 46920486 : z[0] = evaltyp(t_VECSMALL) | evallg(l);
94 46920544 : setlg(x, ly);
95 : #ifdef DEBUG
96 : { long i; for (i = 1; i < l; i++) z[i] = 0; }
97 : #endif
98 : }
99 89140999 : }
100 : /* update lg(z) before affrr(y, z) [ to cater for precision loss ]*/
101 : INLINE void
102 42232967 : affrr_fixlg(GEN y, GEN z) { fixlg(z, lg(y)); affrr(y, z); }
103 :
104 : /*******************************************************************/
105 : /* */
106 : /* ALLOCATE ON STACK */
107 : /* */
108 : /*******************************************************************/
109 : INLINE ulong
110 0 : get_avma(void) { return avma; }
111 : INLINE void
112 >10552*10^7 : set_avma(ulong av) { avma = av; }
113 :
114 : INLINE double
115 300860402 : gc_double(pari_sp av, double d) { set_avma(av); return d; }
116 : INLINE long
117 217182239 : gc_long(pari_sp av, long s) { set_avma(av); return s; }
118 : INLINE ulong
119 26075198 : gc_ulong(pari_sp av, ulong s) { set_avma(av); return s; }
120 : INLINE int
121 43449971 : gc_bool(pari_sp av, int s) { set_avma(av); return s; }
122 : INLINE int
123 1741639 : gc_int(pari_sp av, int s) { set_avma(av); return s; }
124 : INLINE GEN
125 6374196 : gc_NULL(pari_sp av) { set_avma(av); return NULL; }
126 : INLINE GEN
127 12574674757 : gc_const(pari_sp av, GEN x) { set_avma(av); return x; }
128 : INLINE GEN
129 150625 : gc_stoi(pari_sp av, long x) { set_avma(av); return stoi(x); }
130 : INLINE GEN
131 422381 : gc_utoi(pari_sp av, ulong x) { set_avma(av); return utoi(x); }
132 : INLINE GEN
133 1122609 : gc_utoipos(pari_sp av, ulong x) { set_avma(av); return utoipos(x); }
134 :
135 : INLINE GEN
136 76066137915 : new_chunk(size_t x) /* x is a number of longs */
137 : {
138 76066137915 : GEN z = ((GEN) avma) - x;
139 : CHECK_CTRLC
140 76066137915 : if (x > (avma-pari_mainstack->bot) / sizeof(long))
141 14 : new_chunk_resize(x);
142 76066137901 : set_avma((pari_sp)z);
143 : #ifdef MEMSTEP
144 : if (DEBUGMEM>1 && pari_mainstack->memused != DISABLE_MEMUSED) {
145 : long d = (long)pari_mainstack->memused - (long)z;
146 : if (labs(d) > 4*MEMSTEP)
147 : {
148 : pari_mainstack->memused = (pari_sp)z;
149 : err_printf("...%4.0lf Mbytes used\n",
150 : (pari_mainstack->top-pari_mainstack->memused)/1048576.);
151 : }
152 : }
153 : #endif
154 76021485966 : return z;
155 : }
156 :
157 : INLINE char *
158 8504422 : stack_malloc(size_t N)
159 : {
160 8504422 : long n = nchar2nlong(N);
161 8504409 : return (char*)new_chunk(n);
162 : }
163 :
164 : INLINE char *
165 79843923 : stack_malloc_align(size_t N, long k)
166 : {
167 79843923 : ulong d = ((ulong)avma) % k, e = ((ulong)N) % k;
168 79843923 : if (d) (void)new_chunk(d/sizeof(long));
169 79843926 : if (e) N += k-e;
170 79843926 : return (char*) new_chunk(nchar2nlong(N));
171 : }
172 :
173 : INLINE char *
174 85710 : stack_calloc(size_t N)
175 : {
176 85710 : char *p = stack_malloc(N);
177 85710 : memset(p, 0, N); return p;
178 : }
179 :
180 : INLINE char *
181 1094 : stack_calloc_align(size_t N, long k)
182 : {
183 1094 : ulong d = ((ulong)avma) % k, e = ((ulong)N) % k;
184 1094 : if (d) (void)new_chunk(d/sizeof(long));
185 1094 : if (e) N += k-e;
186 1094 : return stack_calloc(N);
187 : }
188 :
189 : /* cgetg(lg(x), typ(x)), set *lx. Implicit unsetisclone() */
190 : INLINE GEN
191 1035172411 : cgetg_copy(GEN x, long *plx) {
192 : GEN y;
193 1035172411 : *plx = lg(x); y = new_chunk((size_t)*plx);
194 1035197627 : y[0] = x[0] & (TYPBITS|LGBITS); return y;
195 : }
196 : INLINE GEN
197 278687 : cgetg_block(long x, long y)
198 : {
199 278687 : GEN z = newblock((size_t)x);
200 278747 : z[0] = CLONEBIT | evaltyp(y) | evallg(x);
201 278701 : return z;
202 : }
203 : INLINE GEN
204 10845986570 : cgetg(long x, long y)
205 : {
206 10845986570 : GEN z = new_chunk((size_t)x);
207 10840232683 : z[0] = evaltyp(y) | evallg(x);
208 10836917587 : return z;
209 : }
210 : INLINE GEN
211 23636339224 : cgeti(long x)
212 : {
213 23636339224 : GEN z = new_chunk((size_t)x);
214 23591800000 : z[0] = evaltyp(t_INT) | evallg(x);
215 23573477453 : return z;
216 : }
217 : INLINE GEN
218 14375175103 : cgetipos(long x)
219 : {
220 14375175103 : GEN z = cgeti(x);
221 14352968047 : z[1] = evalsigne(1) | evallgefint(x);
222 14352968047 : return z;
223 : }
224 : INLINE GEN
225 259522853 : cgetineg(long x)
226 : {
227 259522853 : GEN z = cgeti(x);
228 259524482 : z[1] = evalsigne(-1) | evallgefint(x);
229 259524482 : return z;
230 : }
231 : INLINE GEN
232 39404 : cgetr_block(long x)
233 : {
234 39404 : GEN z = newblock((size_t)x);
235 39429 : z[0] = CLONEBIT | evaltyp(t_REAL) | evallg(x);
236 39424 : return z;
237 : }
238 : INLINE GEN
239 10190346589 : cgetr(long x)
240 : {
241 10190346589 : GEN z = new_chunk((size_t)x);
242 10170281166 : z[0] = evaltyp(t_REAL) | evallg(x);
243 10165040362 : return z;
244 : }
245 :
246 : /*******************************************************************/
247 : /* */
248 : /* COPY, NEGATION, ABSOLUTE VALUE */
249 : /* */
250 : /*******************************************************************/
251 : /* cannot do memcpy because sometimes x and y overlap */
252 : INLINE GEN
253 3117127942 : leafcopy(GEN x)
254 : {
255 3117127942 : long lx = lg(x);
256 3117127942 : GEN y = new_chunk(lx); /* can't use cgetg_copy, in case x,y overlap */
257 17590050707 : while (--lx > 0) y[lx] = x[lx];
258 3116874844 : y[0] = x[0] & (TYPBITS|LGBITS); return y;
259 : }
260 : INLINE GEN
261 7249765141 : icopy(GEN x)
262 : {
263 7249765141 : long i = lgefint(x), lx = i;
264 7249765141 : GEN y = new_chunk(lx); /* can't use cgeti, in case x,y overlap */
265 32180461805 : while (--i > 0) y[i] = x[i];
266 7242482472 : y[0] = evaltyp(t_INT) | evallg(lx);
267 7247612652 : return y;
268 : }
269 : INLINE GEN
270 105938296 : icopyspec(GEN x, long nx)
271 : {
272 105938296 : long i = nx+2, lx = i;
273 105938296 : GEN y = new_chunk(lx); /* can't use cgeti, in case x,y overlap */
274 2819698292 : x -= 2; while (--i >= 2) y[i] = x[i];
275 105937390 : y[1] = evalsigne(1) | evallgefint(lx);
276 105937390 : y[0] = evaltyp(t_INT) | evallg(lx);
277 105937197 : return y;
278 : }
279 770931674 : INLINE GEN rcopy(GEN x) { return leafcopy(x); }
280 707 : INLINE GEN mpcopy(GEN x) { return leafcopy(x); }
281 :
282 : INLINE GEN
283 676147847 : mpabs(GEN x) { GEN y = leafcopy(x); setabssign(y); return y; }
284 : INLINE GEN
285 2612250 : mpabs_shallow(GEN x) { return signe(x) < 0? mpabs(x): x; }
286 643334562 : INLINE GEN absi(GEN x) { return mpabs(x); }
287 53620852 : INLINE GEN absi_shallow(GEN x) { return signe(x) < 0? negi(x): x; }
288 39809 : INLINE GEN absr(GEN x) { return mpabs(x); }
289 :
290 : INLINE GEN
291 795103881 : mpneg(GEN x) { GEN y = leafcopy(x); togglesign(y); return y; }
292 547122029 : INLINE GEN negi(GEN x) { return mpneg(x); }
293 2384672 : INLINE GEN negr(GEN x) { return mpneg(x); }
294 :
295 : /* negate in place */
296 : INLINE void
297 1658836663 : togglesign(GEN x) { if (x[1] & SIGNBITS) { x[1] ^= HIGHBIT; } }
298 : INLINE void
299 740794625 : setabssign(GEN x) { x[1] &= ~HIGHBIT; }
300 : /* negate in place, except universal constants */
301 : INLINE void
302 113950157 : togglesign_safe(GEN *px)
303 : {
304 113950157 : switch(*px - gen_1) /* gen_1, gen_2, gen_m1, gen_m2 */
305 : {
306 2421027 : case 0: *px = gen_m1; break;
307 4 : case 3: *px = gen_m2; break;
308 543921 : case 6: *px = gen_1; break;
309 0 : case 9: *px = gen_2; break;
310 110985205 : default: togglesign(*px);
311 : }
312 113952358 : }
313 : /* setsigne(y, signe(x)) */
314 : INLINE void
315 0 : affectsign(GEN x, GEN y)
316 : {
317 0 : y[1] = (x[1] & SIGNBITS) | (y[1] & ~SIGNBITS);
318 0 : }
319 : /* copies sign in place, except for universal constants */
320 : INLINE void
321 10083561 : affectsign_safe(GEN x, GEN *py)
322 : {
323 10083561 : if (((*py)[1] ^ x[1]) & HIGHBIT) togglesign_safe(py);
324 10083561 : }
325 : /*******************************************************************/
326 : /* */
327 : /* GEN -> LONG, LONG -> GEN */
328 : /* */
329 : /*******************************************************************/
330 : /* assume x != 0, return -x as a t_INT */
331 : INLINE GEN
332 256106345 : utoineg(ulong x) { GEN y = cgetineg(3); y[2] = x; return y; }
333 : /* assume x != 0, return utoi(x) */
334 : INLINE GEN
335 12477125591 : utoipos(ulong x) { GEN y = cgetipos(3); y[2] = x; return y; }
336 : INLINE GEN
337 10600986888 : utoi(ulong x) { return x? utoipos(x): gen_0; }
338 : INLINE GEN
339 722628243 : stoi(long x)
340 : {
341 722628243 : if (!x) return gen_0;
342 508462762 : return x > 0? utoipos((ulong)x): utoineg((ulong)-x);
343 : }
344 :
345 : /* x 2^BIL + y */
346 : INLINE GEN
347 7690410503 : uutoi(ulong x, ulong y)
348 : {
349 : GEN z;
350 7690410503 : if (!x) return utoi(y);
351 790515328 : z = cgetipos(4);
352 793219737 : *int_W_lg(z, 1, 4) = x;
353 793219737 : *int_W_lg(z, 0, 4) = y; return z;
354 : }
355 : /* - (x 2^BIL + y) */
356 : INLINE GEN
357 244617 : uutoineg(ulong x, ulong y)
358 : {
359 : GEN z;
360 244617 : if (!x) return y? utoineg(y): gen_0;
361 10443 : z = cgetineg(4);
362 10443 : *int_W_lg(z, 1, 4) = x;
363 10443 : *int_W_lg(z, 0, 4) = y; return z;
364 : }
365 :
366 : INLINE long
367 484442631 : itos(GEN x)
368 : {
369 484442631 : long s = signe(x);
370 : long u;
371 :
372 484442631 : if (!s) return 0;
373 457618784 : u = x[2];
374 457618784 : if (lgefint(x) > 3 || u < 0)
375 27 : pari_err_OVERFLOW("t_INT-->long assignment");
376 457621403 : return (s>0) ? u : -u;
377 : }
378 : /* as itos, but return 0 if too large. Cf is_bigint */
379 : INLINE long
380 17410057 : itos_or_0(GEN x) {
381 : long n;
382 17410057 : if (lgefint(x) != 3 || (n = x[2]) & HIGHBIT) return 0;
383 15734368 : return signe(x) > 0? n: -n;
384 : }
385 : INLINE ulong
386 158509440 : itou(GEN x)
387 : {
388 158509440 : switch(lgefint(x)) {
389 10818367 : case 2: return 0;
390 147691193 : case 3: return x[2];
391 0 : default:
392 0 : pari_err_OVERFLOW("t_INT-->ulong assignment");
393 : return 0; /* LCOV_EXCL_LINE */
394 : }
395 : }
396 :
397 : /* as itou, but return 0 if too large. Cf is_bigint */
398 : INLINE ulong
399 2755226 : itou_or_0(GEN x) {
400 2755226 : if (lgefint(x) != 3) return 0;
401 2742296 : return (ulong)x[2];
402 : }
403 :
404 : INLINE ulong
405 5477034 : umuluu_or_0(ulong x, ulong y)
406 : {
407 : ulong z;
408 : LOCAL_HIREMAINDER;
409 5477034 : z = mulll(x, y);
410 5477034 : return hiremainder? 0: z;
411 : }
412 : /* return x*y if <= n, else 0. Beware overflow */
413 : INLINE ulong
414 5645062 : umuluu_le(ulong x, ulong y, ulong n)
415 : {
416 : ulong z;
417 : LOCAL_HIREMAINDER;
418 5645062 : z = mulll(x, y);
419 5645062 : return (hiremainder || z > n)? 0: z;
420 : }
421 :
422 : INLINE GEN
423 108258214 : real_0_bit(long bitprec) { GEN x=cgetr(2); x[1]=evalexpo(bitprec); return x; }
424 : INLINE GEN
425 549287 : real_0(long prec) { return real_0_bit(-prec2nbits(prec)); }
426 : INLINE GEN
427 2711272 : real_1_bit(long bit) { return real_1(nbits2prec(bit)); }
428 : INLINE GEN
429 87158746 : real_1(long prec) {
430 87158746 : GEN x = cgetr(prec);
431 : long i;
432 87150212 : x[1] = evalsigne(1) | _evalexpo(0);
433 446307404 : x[2] = (long)HIGHBIT; for (i=3; i<prec; i++) x[i] = 0;
434 87150212 : return x;
435 : }
436 : INLINE GEN
437 322 : real_m1(long prec) {
438 322 : GEN x = cgetr(prec);
439 : long i;
440 322 : x[1] = evalsigne(-1) | _evalexpo(0);
441 1601 : x[2] = (long)HIGHBIT; for (i=3; i<prec; i++) x[i] = 0;
442 322 : return x;
443 : }
444 :
445 : /* 2.^n */
446 : INLINE GEN
447 602941 : real2n(long n, long prec) { GEN z = real_1(prec); setexpo(z, n); return z; }
448 : INLINE GEN
449 0 : real_m2n(long n, long prec) { GEN z = real_m1(prec); setexpo(z, n); return z; }
450 : INLINE GEN
451 419609607 : stor(long s, long prec) { GEN z = cgetr(prec); affsr(s,z); return z; }
452 : INLINE GEN
453 12410006 : utor(ulong s, long prec){ GEN z = cgetr(prec); affur(s,z); return z; }
454 : INLINE GEN
455 703867093 : itor(GEN x, long prec) { GEN z = cgetr(prec); affir(x,z); return z; }
456 : INLINE GEN
457 229552182 : rtor(GEN x, long prec) { GEN z = cgetr(prec); affrr(x,z); return z; }
458 :
459 16871434 : INLINE ulong int_bit(GEN x, long n)
460 : {
461 16871434 : long r, q = dvmdsBIL(n, &r);
462 16867058 : return q < lgefint(x)-2?((ulong)*int_W(x,q) >> r) & 1UL:0;
463 : }
464 :
465 : /*******************************************************************/
466 : /* */
467 : /* COMPARISON */
468 : /* */
469 : /*******************************************************************/
470 : INLINE int
471 1195345 : cmpss(long a, long b)
472 1195345 : { return a>b? 1: (a<b? -1: 0); }
473 :
474 : INLINE int
475 1401863031 : cmpuu(ulong a, ulong b)
476 1401863031 : { return a>b? 1: (a<b? -1: 0); }
477 :
478 : INLINE int
479 1690787 : cmpir(GEN x, GEN y)
480 : {
481 : pari_sp av;
482 : GEN z;
483 :
484 1690787 : if (!signe(x)) return -signe(y);
485 475658 : if (!signe(y))
486 : {
487 2241 : if (expo(y) >= expi(x)) return 0;
488 2213 : return signe(x);
489 : }
490 473417 : av=avma; z = itor(x, realprec(y)); set_avma(av);
491 473408 : return cmprr(z,y); /* cmprr does no memory adjustment */
492 : }
493 : INLINE int
494 409272 : cmpri(GEN x, GEN y) { return -cmpir(y,x); }
495 : INLINE int
496 606130 : cmpsr(long x, GEN y)
497 : {
498 : pari_sp av;
499 : GEN z;
500 :
501 606130 : if (!x) return -signe(y);
502 606130 : av=avma; z = stor(x, LOWDEFAULTPREC); set_avma(av);
503 606130 : return cmprr(z,y);
504 : }
505 : INLINE int
506 40989 : cmprs(GEN x, long y) { return -cmpsr(y,x); }
507 : /* compare x and y */
508 : INLINE int
509 8636550 : cmpui(ulong x, GEN y)
510 : {
511 : ulong p;
512 8636550 : if (!x) return -signe(y);
513 8636487 : if (signe(y) <= 0) return 1;
514 8524718 : if (lgefint(y) > 3) return -1;
515 8461573 : p = y[2]; if (p == x) return 0;
516 8392963 : return p < x ? 1 : -1;
517 : }
518 : INLINE int
519 8632381 : cmpiu(GEN x, ulong y) { return -cmpui(y,x); }
520 : /* compare x and |y| */
521 : INLINE int
522 33797515 : abscmpui(ulong x, GEN y)
523 : {
524 33797515 : long l = lgefint(y);
525 : ulong p;
526 :
527 33797515 : if (!x) return (l > 2)? -1: 0;
528 33797501 : if (l == 2) return 1;
529 33423804 : if (l > 3) return -1;
530 33396979 : p = y[2]; if (p == x) return 0;
531 32872943 : return p < x ? 1 : -1;
532 : }
533 : INLINE int
534 33796051 : abscmpiu(GEN x, ulong y) { return -abscmpui(y,x); }
535 : INLINE int
536 3852492 : cmpsi(long x, GEN y)
537 : {
538 : ulong p;
539 :
540 3852492 : if (!x) return -signe(y);
541 :
542 3851156 : if (x > 0)
543 : {
544 3850176 : if (signe(y)<=0) return 1;
545 3849882 : if (lgefint(y)>3) return -1;
546 3832083 : p = y[2]; if (p == (ulong)x) return 0;
547 3766731 : return p < (ulong)x ? 1 : -1;
548 : }
549 :
550 980 : if (signe(y)>=0) return -1;
551 119 : if (lgefint(y)>3) return 1;
552 119 : p = y[2]; if (p == (ulong)-x) return 0;
553 14 : return p < (ulong)(-x) ? -1 : 1;
554 : }
555 : INLINE int
556 3830611 : cmpis(GEN x, long y) { return -cmpsi(y,x); }
557 : INLINE int
558 2109072 : mpcmp(GEN x, GEN y)
559 : {
560 2109072 : if (typ(x)==t_INT)
561 60362 : return (typ(y)==t_INT) ? cmpii(x,y) : cmpir(x,y);
562 2048710 : return (typ(y)==t_INT) ? -cmpir(y,x) : cmprr(x,y);
563 : }
564 :
565 : /* x == y ? */
566 : INLINE int
567 2891746 : equalui(ulong x, GEN y)
568 : {
569 2891746 : if (!x) return !signe(y);
570 2891053 : if (signe(y) <= 0 || lgefint(y) != 3) return 0;
571 2880150 : return ((ulong)y[2] == (ulong)x);
572 : }
573 : /* x == y ? */
574 : INLINE int
575 488728 : equalsi(long x, GEN y)
576 : {
577 488728 : if (!x) return !signe(y);
578 488728 : if (x > 0)
579 : {
580 486530 : if (signe(y) <= 0 || lgefint(y) != 3) return 0;
581 434981 : return ((ulong)y[2] == (ulong)x);
582 : }
583 2198 : if (signe(y) >= 0 || lgefint(y) != 3) return 0;
584 2078 : return ((ulong)y[2] == (ulong)-x);
585 : }
586 : /* x == |y| ? */
587 : INLINE int
588 42113888 : absequalui(ulong x, GEN y)
589 : {
590 42113888 : if (!x) return !signe(y);
591 42113888 : return (lgefint(y) == 3 && (ulong)y[2] == x);
592 : }
593 : INLINE int
594 40423458 : absequaliu(GEN x, ulong y) { return absequalui(y,x); }
595 : INLINE int
596 488546 : equalis(GEN x, long y) { return equalsi(y,x); }
597 : INLINE int
598 2891742 : equaliu(GEN x, ulong y) { return equalui(y,x); }
599 :
600 : /* assume x != 0, is |x| == 2^n ? */
601 : INLINE int
602 1149528 : absrnz_equal2n(GEN x) {
603 1149528 : if ((ulong)x[2]==HIGHBIT)
604 : {
605 32575 : long i, lx = lg(x);
606 105038 : for (i = 3; i < lx; i++)
607 79044 : if (x[i]) return 0;
608 25994 : return 1;
609 : }
610 1116953 : return 0;
611 : }
612 : /* assume x != 0, is |x| == 1 ? */
613 : INLINE int
614 3658397 : absrnz_equal1(GEN x) { return !expo(x) && absrnz_equal2n(x); }
615 :
616 : INLINE long
617 8604531245 : maxss(long x, long y) { return x>y?x:y; }
618 : INLINE long
619 1525383386 : minss(long x, long y) { return x<y?x:y; }
620 : INLINE long
621 5189888 : minuu(ulong x, ulong y) { return x<y?x:y; }
622 : INLINE long
623 9708811 : maxuu(ulong x, ulong y) { return x>y?x:y; }
624 : INLINE double
625 2796252 : maxdd(double x, double y) { return x>y?x:y; }
626 : INLINE double
627 592757 : mindd(double x, double y) { return x<y?x:y; }
628 :
629 : /*******************************************************************/
630 : /* */
631 : /* ADD / SUB */
632 : /* */
633 : /*******************************************************************/
634 : INLINE GEN
635 25067 : subuu(ulong x, ulong y)
636 : {
637 : ulong z;
638 : LOCAL_OVERFLOW;
639 25067 : z = subll(x, y);
640 25067 : return overflow? utoineg(-z): utoi(z);
641 : }
642 : INLINE GEN
643 3175987798 : adduu(ulong x, ulong y) { ulong t = x+y; return uutoi((t < x), t); }
644 :
645 : INLINE GEN
646 25067 : addss(long x, long y)
647 : {
648 25067 : if (!x) return stoi(y);
649 25067 : if (!y) return stoi(x);
650 25067 : if (x > 0) return y > 0? adduu(x,y): subuu(x, -y);
651 :
652 25067 : if (y > 0) return subuu(y, -x);
653 : else { /* - adduu(-x, -y) */
654 0 : ulong t = (-x)+(-y); return uutoineg((t < (ulong)(-x)), t);
655 : }
656 : }
657 25067 : INLINE GEN subss(long x, long y) { return addss(-y,x); }
658 :
659 : INLINE GEN
660 7590668763 : subii(GEN x, GEN y)
661 : {
662 7590668763 : if (x==y) return gen_0; /* frequent with x = y = gen_0 */
663 5855043712 : return addii_sign(x, signe(x), y, -signe(y));
664 : }
665 : INLINE GEN
666 9532340291 : addii(GEN x, GEN y) { return addii_sign(x, signe(x), y, signe(y)); }
667 : INLINE GEN
668 2091269237 : addrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, signe(y)); }
669 : INLINE GEN
670 758927227 : subrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, -signe(y)); }
671 : INLINE GEN
672 449321182 : addir(GEN x, GEN y) { return addir_sign(x, signe(x), y, signe(y)); }
673 : INLINE GEN
674 2807663 : subir(GEN x, GEN y) { return addir_sign(x, signe(x), y, -signe(y)); }
675 : INLINE GEN
676 7168524 : subri(GEN x, GEN y) { return addir_sign(y, -signe(y), x, signe(x)); }
677 : INLINE GEN
678 284559821 : addsi(long x, GEN y) { return addsi_sign(x, y, signe(y)); }
679 : INLINE GEN
680 94466855 : addui(ulong x, GEN y) { return addui_sign(x, y, signe(y)); }
681 : INLINE GEN
682 5609758 : subsi(long x, GEN y) { return addsi_sign(x, y, -signe(y)); }
683 : INLINE GEN
684 125309618 : subui(ulong x, GEN y) { return addui_sign(x, y, -signe(y)); }
685 :
686 : /*******************************************************************/
687 : /* */
688 : /* MOD, REM, DIV */
689 : /* */
690 : /*******************************************************************/
691 87870044 : INLINE ulong mod2BIL(GEN x) { return *int_LSW(x); }
692 0 : INLINE long mod64(GEN x) { return mod2BIL(x) & 63; }
693 259 : INLINE long mod32(GEN x) { return mod2BIL(x) & 31; }
694 235247 : INLINE long mod16(GEN x) { return mod2BIL(x) & 15; }
695 12422584 : INLINE long mod8(GEN x) { return mod2BIL(x) & 7; }
696 4076827 : INLINE long mod4(GEN x) { return mod2BIL(x) & 3; }
697 51997547 : INLINE long mod2(GEN x) { return mod2BIL(x) & 1; }
698 : INLINE int
699 83339100 : mpodd(GEN x) { return signe(x) && mod2(x); }
700 : /* x mod 2^n, n < BITS_IN_LONG */
701 : INLINE ulong
702 20161037 : umodi2n(GEN x, long n)
703 : {
704 20161037 : long s = signe(x);
705 20161037 : const ulong _2n = 1UL << n;
706 : ulong m;
707 20161037 : if (!s) return 0;
708 20146302 : m = *int_LSW(x) & (_2n - 1);
709 20146302 : if (s < 0 && m) m = _2n - m;
710 20146302 : return m;
711 : }
712 0 : INLINE ulong Mod64(GEN x){ return umodi2n(x,6); }
713 167642 : INLINE ulong Mod32(GEN x){ return umodi2n(x,5); }
714 244224 : INLINE ulong Mod16(GEN x){ return umodi2n(x,4); }
715 1780919 : INLINE ulong Mod8(GEN x) { return umodi2n(x,3); }
716 16278317 : INLINE ulong Mod4(GEN x) { return umodi2n(x,2); }
717 1689190 : INLINE ulong Mod2(GEN x) { return umodi2n(x,1); }
718 :
719 : INLINE GEN
720 56193835 : truedivii(GEN a,GEN b) { return truedvmdii(a,b,NULL); }
721 : INLINE GEN
722 237085 : truedivis(GEN a, long b) { return truedvmdis(a,b,NULL); }
723 : INLINE GEN
724 6197775 : truedivsi(long a, GEN b) { return truedvmdsi(a,b,NULL); }
725 :
726 : INLINE GEN
727 12490104 : divii(GEN a, GEN b) { return dvmdii(a,b,NULL); }
728 : INLINE GEN
729 2444972927 : remii(GEN a, GEN b) { return dvmdii(a,b,ONLY_REM); }
730 :
731 : INLINE GEN
732 0 : divss(long x, long y) { return stoi(x / y); }
733 : INLINE GEN
734 0 : modss(long x, long y) { return utoi(smodss(x, y)); }
735 : INLINE GEN
736 0 : remss(long x, long y) { return stoi(x % y); }
737 : INLINE long
738 6418 : smodss(long x, long y)
739 : {
740 6418 : long r = x%y;
741 6418 : return (r >= 0)? r: labs(y) + r;
742 : }
743 : INLINE ulong
744 713939470 : umodsu(long x, ulong y)
745 : {
746 713939470 : return x>=0 ? x%y: Fl_neg((-x)%y, y);
747 : }
748 :
749 : INLINE long
750 0 : sdivss_rem(long x, long y, long *r)
751 : {
752 : long q;
753 : LOCAL_HIREMAINDER;
754 0 : if (!y) pari_err_INV("sdivss_rem",gen_0);
755 0 : hiremainder = 0; q = divll((ulong)labs(x),(ulong)labs(y));
756 0 : if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
757 0 : if (y < 0) q = -q;
758 0 : *r = hiremainder; return q;
759 : }
760 : INLINE GEN
761 0 : divss_rem(long x, long y, long *r) { return stoi(sdivss_rem(x,y,r)); }
762 : INLINE ulong
763 158379529 : udivuu_rem(ulong x, ulong y, ulong *r)
764 : {
765 158379529 : if (!y) pari_err_INV("udivuu_rem",gen_0);
766 158379529 : *r = x % y; return x / y;
767 : }
768 : INLINE ulong
769 1870395 : ceildivuu(ulong a, ulong b)
770 : {
771 : ulong c;
772 1870395 : if (!a) return 0;
773 1692981 : c = a / b; return (a % b)? c+1: c;
774 : }
775 :
776 : INLINE ulong
777 13295 : uabsdivui_rem(ulong x, GEN y, ulong *r)
778 : {
779 13295 : long q, s = signe(y);
780 : LOCAL_HIREMAINDER;
781 :
782 13295 : if (!s) pari_err_INV("uabsdivui_rem",gen_0);
783 13295 : if (!x || lgefint(y)>3) { *r = x; return 0; }
784 12634 : hiremainder=0; q = (long)divll(x, (ulong)y[2]);
785 12634 : if (s < 0) q = -q;
786 12634 : *r = hiremainder; return q;
787 : }
788 :
789 : /* assume d != 0 and |n| / d can be represented as an ulong.
790 : * Return |n|/d, set *r = |n| % d */
791 : INLINE ulong
792 10187366 : uabsdiviu_rem(GEN n, ulong d, ulong *r)
793 : {
794 10187366 : switch(lgefint(n))
795 : {
796 0 : case 2: *r = 0; return 0;
797 10187366 : case 3:
798 : {
799 10187366 : ulong nn = n[2];
800 10187366 : *r = nn % d; return nn / d;
801 : }
802 0 : default: /* 4 */
803 : {
804 : ulong n1, n0, q;
805 : LOCAL_HIREMAINDER;
806 0 : n0 = *int_W(n,0);
807 0 : n1 = *int_W(n,1);
808 0 : hiremainder = n1;
809 0 : q = divll(n0, d);
810 0 : *r = hiremainder; return q;
811 : }
812 : }
813 : }
814 :
815 : INLINE long
816 51288571 : sdivsi_rem(long x, GEN y, long *r)
817 : {
818 51288571 : long q, s = signe(y);
819 : LOCAL_HIREMAINDER;
820 :
821 51288571 : if (!s) pari_err_INV("sdivsi_rem",gen_0);
822 51288572 : if (!x || lgefint(y)>3 || ((long)y[2]) < 0) { *r = x; return 0; }
823 49334326 : hiremainder=0; q = (long)divll(labs(x), (ulong)y[2]);
824 49334326 : if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
825 49334326 : if (s < 0) q = -q;
826 49334326 : *r = hiremainder; return q;
827 : }
828 : INLINE GEN
829 0 : divsi_rem(long s, GEN y, long *r) { return stoi(sdivsi_rem(s,y,r)); }
830 :
831 : INLINE long
832 100946 : sdivsi(long x, GEN y)
833 : {
834 100946 : long q, s = signe(y);
835 :
836 100946 : if (!s) pari_err_INV("sdivsi",gen_0);
837 100946 : if (!x || lgefint(y)>3 || ((long)y[2]) < 0) return 0;
838 100841 : q = labs(x) / y[2];
839 100841 : if (x < 0) q = -q;
840 100841 : if (s < 0) q = -q;
841 100841 : return q;
842 : }
843 :
844 : INLINE GEN
845 0 : dvmdss(long x, long y, GEN *z)
846 : {
847 : long r;
848 0 : GEN q = divss_rem(x,y, &r);
849 0 : *z = stoi(r); return q;
850 : }
851 : INLINE long
852 5501389687 : dvmdsBIL(long n, long *r) { *r = remsBIL(n); return divsBIL(n); }
853 : INLINE ulong
854 165473043 : dvmduBIL(ulong n, ulong *r) { *r = remsBIL(n); return divsBIL(n); }
855 : INLINE GEN
856 0 : dvmdsi(long x, GEN y, GEN *z)
857 : {
858 : long r;
859 0 : GEN q = divsi_rem(x,y, &r);
860 0 : *z = stoi(r); return q;
861 : }
862 : INLINE GEN
863 0 : dvmdis(GEN x, long y, GEN *z)
864 : {
865 : long r;
866 0 : GEN q = divis_rem(x,y, &r);
867 0 : *z = stoi(r); return q;
868 : }
869 :
870 : INLINE long
871 21005043 : smodis(GEN x, long y)
872 : {
873 21005043 : pari_sp av = avma;
874 21005043 : long r; (void)divis_rem(x,y, &r);
875 21005043 : return gc_long(av, (r >= 0)? r: labs(y) + r);
876 : }
877 : INLINE GEN
878 19602337 : modis(GEN x, long y) { return stoi(smodis(x,y)); }
879 : INLINE GEN
880 45090552 : modsi(long x, GEN y) {
881 45090552 : long r; (void)sdivsi_rem(x, y, &r);
882 45090552 : return (r >= 0)? stoi(r): addsi_sign(r, y, 1);
883 : }
884 :
885 : INLINE ulong
886 1339528 : umodui(ulong x, GEN y)
887 : {
888 1339528 : if (!signe(y)) pari_err_INV("umodui",gen_0);
889 1339528 : if (!x || lgefint(y) > 3) return x;
890 415201 : return x % (ulong)y[2];
891 : }
892 :
893 : INLINE ulong
894 210600 : ugcdiu(GEN x, ulong y) { return ugcd(umodiu(x,y), y); }
895 : INLINE ulong
896 2737 : ugcdui(ulong y, GEN x) { return ugcd(umodiu(x,y), y); }
897 :
898 : INLINE GEN
899 0 : remsi(long x, GEN y)
900 0 : { long r; (void)sdivsi_rem(x,y, &r); return stoi(r); }
901 : INLINE GEN
902 0 : remis(GEN x, long y)
903 : {
904 0 : pari_sp av = avma;
905 : long r;
906 0 : (void)divis_rem(x,y, &r); return gc_stoi(av, r);
907 : }
908 :
909 : INLINE GEN
910 0 : rdivis(GEN x, long y, long prec)
911 : {
912 0 : GEN z = cgetr(prec);
913 0 : pari_sp av = avma;
914 0 : affrr(divrs(itor(x,prec), y),z);
915 0 : set_avma(av); return z;
916 : }
917 : INLINE GEN
918 0 : rdivsi(long x, GEN y, long prec)
919 : {
920 0 : GEN z = cgetr(prec);
921 0 : pari_sp av = avma;
922 0 : affrr(divsr(x, itor(y,prec)), z);
923 0 : set_avma(av); return z;
924 : }
925 : INLINE GEN
926 839647 : rdivss(long x, long y, long prec)
927 : {
928 839647 : GEN z = cgetr(prec);
929 839647 : pari_sp av = avma;
930 839647 : affrr(divrs(stor(x, prec), y), z);
931 839647 : set_avma(av); return z;
932 : }
933 :
934 : INLINE void
935 7803971 : rdiviiz(GEN x, GEN y, GEN z)
936 : {
937 7803971 : long prec = realprec(z), lx = lgefint(x), ly = lgefint(y);
938 7803971 : if (lx == 2) { affur(0, z); return; }
939 7803971 : if (ly == 3)
940 : {
941 2215242 : affir(x, z); if (signe(y) < 0) togglesign(z);
942 2215190 : affrr(divru(z, y[2]), z);
943 : }
944 5588729 : else if (lx > prec + 1 || ly > prec + 1)
945 : {
946 1851419 : affir(x,z); affrr(divri(z, y), z);
947 : }
948 : else
949 : {
950 3737310 : long b = bit_accuracy(prec) + expi(y) - expi(x) + 1;
951 3737397 : GEN q = divii(b > 0? shifti(x, b): x, y);
952 3737372 : affir(q, z); if (b > 0) shiftr_inplace(z, -b);
953 : }
954 7804085 : set_avma((ulong)z);
955 : }
956 : INLINE GEN
957 7764354 : rdivii(GEN x, GEN y, long prec)
958 7764354 : { GEN z = cgetr(prec); rdiviiz(x, y, z); return z; }
959 : INLINE GEN
960 7350733 : fractor(GEN x, long prec)
961 7350733 : { return rdivii(gel(x,1), gel(x,2), prec); }
962 :
963 : INLINE int
964 16044756 : dvdii(GEN x, GEN y)
965 : {
966 16044756 : pari_sp av = avma;
967 : GEN r;
968 16044756 : if (!signe(x)) return 1;
969 14756837 : if (!signe(y)) return 0;
970 14756837 : r = remii(x,y);
971 14763723 : return gc_bool(av, r == gen_0);
972 : }
973 : INLINE int
974 371 : dvdsi(long x, GEN y)
975 : {
976 371 : if (x == 0) return 1;
977 266 : if (!signe(y)) return 0;
978 266 : if (lgefint(y) != 3) return 0;
979 259 : return x % y[2] == 0;
980 : }
981 : INLINE int
982 167181 : dvdui(ulong x, GEN y)
983 : {
984 167181 : if (x == 0) return 1;
985 167181 : if (!signe(y)) return 0;
986 167181 : if (lgefint(y) != 3) return 0;
987 156574 : return x % y[2] == 0;
988 : }
989 : INLINE int
990 33331 : dvdis(GEN x, long y)
991 33331 : { return y? smodis(x, y) == 0: signe(x) == 0; }
992 : INLINE int
993 573605 : dvdiu(GEN x, ulong y)
994 573605 : { return y? umodiu(x, y) == 0: signe(x) == 0; }
995 :
996 : INLINE int
997 0 : dvdisz(GEN x, long y, GEN z)
998 : {
999 0 : const pari_sp av = avma;
1000 : long r;
1001 0 : GEN p1 = divis_rem(x,y, &r);
1002 0 : set_avma(av); if (r) return 0;
1003 0 : affii(p1,z); return 1;
1004 : }
1005 : INLINE int
1006 0 : dvdiuz(GEN x, ulong y, GEN z)
1007 : {
1008 0 : const pari_sp av = avma;
1009 : ulong r;
1010 0 : GEN p1 = absdiviu_rem(x,y, &r);
1011 0 : set_avma(av); if (r) return 0;
1012 0 : affii(p1,z); return 1;
1013 : }
1014 : INLINE int
1015 6097 : dvdiiz(GEN x, GEN y, GEN z)
1016 : {
1017 6097 : const pari_sp av=avma;
1018 6097 : GEN p2, p1 = dvmdii(x,y,&p2);
1019 6097 : if (signe(p2)) return gc_bool(av,0);
1020 4896 : affii(p1,z); return gc_bool(av,1);
1021 : }
1022 :
1023 : INLINE ulong
1024 84375408 : remlll_pre(ulong u2, ulong u1, ulong u0, ulong n, ulong ninv)
1025 : {
1026 84375408 : u1 = remll_pre(u2, u1, n, ninv);
1027 85371466 : return remll_pre(u1, u0, n, ninv);
1028 : }
1029 :
1030 : INLINE ulong
1031 1834711500 : Fl_sqr_pre(ulong a, ulong p, ulong pi)
1032 : {
1033 : ulong x;
1034 : LOCAL_HIREMAINDER;
1035 1834711500 : x = mulll(a,a);
1036 1834711500 : return remll_pre(hiremainder, x, p, pi);
1037 : }
1038 :
1039 : INLINE ulong
1040 3142425219 : Fl_mul_pre(ulong a, ulong b, ulong p, ulong pi)
1041 : {
1042 : ulong x;
1043 : LOCAL_HIREMAINDER;
1044 3142425219 : x = mulll(a,b);
1045 3142425219 : return remll_pre(hiremainder, x, p, pi);
1046 : }
1047 :
1048 : INLINE ulong
1049 6424350715 : Fl_addmul_pre(ulong y0, ulong x0, ulong x1, ulong p, ulong pi)
1050 : {
1051 : ulong l0, h0;
1052 : LOCAL_HIREMAINDER;
1053 6424350715 : hiremainder = y0;
1054 6424350715 : l0 = addmul(x0, x1); h0 = hiremainder;
1055 6424350715 : return remll_pre(h0, l0, p, pi);
1056 : }
1057 :
1058 : INLINE ulong
1059 53603709 : Fl_addmulmul_pre(ulong x0, ulong y0, ulong x1, ulong y1, ulong p, ulong pi)
1060 : {
1061 : ulong l0, l1, h0, h1;
1062 : LOCAL_OVERFLOW;
1063 : LOCAL_HIREMAINDER;
1064 53603709 : l0 = mulll(x0, y0); h0 = hiremainder;
1065 53603709 : l1 = mulll(x1, y1); h1 = hiremainder;
1066 53603709 : l0 = addll(l0, l1); h0 = addllx(h0, h1);
1067 53603709 : return overflow ? remlll_pre(1, h0, l0, p, pi): remll_pre(h0, l0, p, pi);
1068 : }
1069 :
1070 : INLINE ulong
1071 199908 : Fl_ellj_pre(ulong a4, ulong a6, ulong p, ulong pi)
1072 : {
1073 : /* a43 = 4 a4^3 */
1074 199908 : ulong a43 = Fl_double(Fl_double(
1075 : Fl_mul_pre(a4, Fl_sqr_pre(a4, p, pi), p, pi), p), p);
1076 : /* a62 = 27 a6^2 */
1077 199912 : ulong a62 = Fl_mul_pre(Fl_sqr_pre(a6, p, pi), 27 % p, p, pi);
1078 199916 : ulong z1 = Fl_mul_pre(a43, 1728 % p, p, pi);
1079 199916 : ulong z2 = Fl_add(a43, a62, p);
1080 199916 : return Fl_div(z1, z2, p);
1081 : }
1082 :
1083 : /*******************************************************************/
1084 : /* */
1085 : /* MP (INT OR REAL) */
1086 : /* */
1087 : /*******************************************************************/
1088 : INLINE GEN
1089 42 : mptrunc(GEN x) { return typ(x)==t_INT? icopy(x): truncr(x); }
1090 : INLINE GEN
1091 0 : mpfloor(GEN x) { return typ(x)==t_INT? icopy(x): floorr(x); }
1092 : INLINE GEN
1093 0 : mpceil(GEN x) { return typ(x)==t_INT? icopy(x): ceilr(x); }
1094 : INLINE GEN
1095 1416579 : mpround(GEN x) { return typ(x) == t_INT? icopy(x): roundr(x); }
1096 :
1097 : INLINE long
1098 636513 : mpexpo(GEN x) { return typ(x) == t_INT? expi(x): expo(x); }
1099 :
1100 : INLINE GEN
1101 97918401 : mpadd(GEN x, GEN y)
1102 : {
1103 97918401 : if (typ(x)==t_INT)
1104 18097270 : return (typ(y)==t_INT) ? addii(x,y) : addir(x,y);
1105 79821131 : return (typ(y)==t_INT) ? addir(y,x) : addrr(x,y);
1106 : }
1107 : INLINE GEN
1108 53330841 : mpsub(GEN x, GEN y)
1109 : {
1110 53330841 : if (typ(x)==t_INT)
1111 506829 : return (typ(y)==t_INT) ? subii(x,y) : subir(x,y);
1112 52824012 : return (typ(y)==t_INT) ? subri(x,y) : subrr(x,y);
1113 : }
1114 : INLINE GEN
1115 153258273 : mpmul(GEN x, GEN y)
1116 : {
1117 153258273 : if (typ(x)==t_INT)
1118 35089732 : return (typ(y)==t_INT) ? mulii(x,y) : mulir(x,y);
1119 118168541 : return (typ(y)==t_INT) ? mulir(y,x) : mulrr(x,y);
1120 : }
1121 : INLINE GEN
1122 18066123 : mpsqr(GEN x) { return (typ(x)==t_INT) ? sqri(x) : sqrr(x); }
1123 : INLINE GEN
1124 624488 : mpdiv(GEN x, GEN y)
1125 : {
1126 624488 : if (typ(x)==t_INT)
1127 235759 : return (typ(y)==t_INT) ? divii(x,y) : divir(x,y);
1128 388729 : return (typ(y)==t_INT) ? divri(x,y) : divrr(x,y);
1129 : }
1130 :
1131 : /*******************************************************************/
1132 : /* */
1133 : /* Z/nZ, n ULONG */
1134 : /* */
1135 : /*******************************************************************/
1136 : INLINE ulong
1137 416393397 : Fl_double(ulong a, ulong p)
1138 : {
1139 416393397 : ulong res = a << 1;
1140 416393397 : return (res >= p || res < a) ? res - p : res;
1141 : }
1142 : INLINE ulong
1143 84353966 : Fl_triple(ulong a, ulong p)
1144 : {
1145 84353966 : ulong res = a << 1;
1146 84353966 : if (res >= p || res < a) res -= p;
1147 84353966 : res += a;
1148 84353966 : return (res >= p || res < a)? res - p: res;
1149 : }
1150 : INLINE ulong
1151 16174714 : Fl_halve(ulong a, ulong p)
1152 : {
1153 : ulong ap, ap2;
1154 16174714 : if ((a&1UL)==0) return a>>1;
1155 8109110 : ap = a + p; ap2 = ap>>1;
1156 8109110 : return ap>=a ? ap2: (ap2|HIGHBIT);
1157 : }
1158 :
1159 : INLINE ulong
1160 3917885655 : Fl_add(ulong a, ulong b, ulong p)
1161 : {
1162 3917885655 : ulong res = a + b;
1163 3917885655 : return (res >= p || res < a) ? res - p : res;
1164 : }
1165 : INLINE ulong
1166 685565803 : Fl_neg(ulong x, ulong p) { return x ? p - x: 0; }
1167 :
1168 : INLINE ulong
1169 6706079859 : Fl_sub(ulong a, ulong b, ulong p)
1170 : {
1171 6706079859 : ulong res = a - b;
1172 6706079859 : return (res > a) ? res + p: res;
1173 : }
1174 :
1175 : /* centerlift(u mod p) */
1176 : INLINE long
1177 3894921 : Fl_center(ulong u, ulong p, ulong ps2) { return (long) (u > ps2)? u - p: u; }
1178 :
1179 : INLINE ulong
1180 2246844973 : Fl_mul(ulong a, ulong b, ulong p)
1181 : {
1182 : ulong x;
1183 : LOCAL_HIREMAINDER;
1184 2246844973 : x = mulll(a,b);
1185 2246844973 : if (!hiremainder) return x % p;
1186 360283686 : (void)divll(x,p); return hiremainder;
1187 : }
1188 : INLINE ulong
1189 92439947 : Fl_sqr(ulong a, ulong p)
1190 : {
1191 : ulong x;
1192 : LOCAL_HIREMAINDER;
1193 92439947 : x = mulll(a,a);
1194 92439947 : if (!hiremainder) return x % p;
1195 25698467 : (void)divll(x,p); return hiremainder;
1196 : }
1197 : /* don't assume that p is prime: can't special case a = 0 */
1198 : INLINE ulong
1199 31918036 : Fl_div(ulong a, ulong b, ulong p)
1200 31918036 : { return Fl_mul(a, Fl_inv(b, p), p); }
1201 :
1202 : /*******************************************************************/
1203 : /* */
1204 : /* DEFINED FROM EXISTING ONE EXPLOITING COMMUTATIVITY */
1205 : /* */
1206 : /*******************************************************************/
1207 : INLINE GEN
1208 1096067 : addri(GEN x, GEN y) { return addir(y,x); }
1209 : INLINE GEN
1210 149980143 : addis(GEN x, long s) { return addsi(s,x); }
1211 : INLINE GEN
1212 91127114 : addiu(GEN x, ulong s) { return addui(s,x); }
1213 : INLINE GEN
1214 10311942 : addrs(GEN x, long s) { return addsr(s,x); }
1215 :
1216 : INLINE GEN
1217 120190441 : subiu(GEN x, long y) { GEN z = subui(y, x); togglesign(z); return z; }
1218 : INLINE GEN
1219 170064 : subis(GEN x, long y) { return addsi(-y,x); }
1220 : INLINE GEN
1221 12743298 : subrs(GEN x, long y) { return addsr(-y,x); }
1222 :
1223 : INLINE GEN
1224 415683579 : mulis(GEN x, long s) { return mulsi(s,x); }
1225 : INLINE GEN
1226 355019853 : muliu(GEN x, ulong s) { return mului(s,x); }
1227 : INLINE GEN
1228 2680835 : mulru(GEN x, ulong s) { return mulur(s,x); }
1229 : INLINE GEN
1230 36008825 : mulri(GEN x, GEN s) { return mulir(s,x); }
1231 : INLINE GEN
1232 10271246 : mulrs(GEN x, long s) { return mulsr(s,x); }
1233 :
1234 : /*******************************************************************/
1235 : /* */
1236 : /* VALUATION, EXPONENT, SHIFTS */
1237 : /* */
1238 : /*******************************************************************/
1239 : INLINE long
1240 164914021 : vali(GEN x)
1241 : {
1242 : long i;
1243 : GEN xp;
1244 :
1245 164914021 : if (!signe(x)) return -1;
1246 164837449 : xp=int_LSW(x);
1247 169749604 : for (i=0; !*xp; i++) xp=int_nextW(xp);
1248 164837449 : return vals(*xp) + i * BITS_IN_LONG;
1249 : }
1250 :
1251 : /* assume x > 0 */
1252 : INLINE long
1253 781390061 : expu(ulong x) { return (BITS_IN_LONG-1) - (long)bfffo(x); }
1254 :
1255 : INLINE long
1256 1418178901 : expi(GEN x)
1257 : {
1258 1418178901 : const long lx=lgefint(x);
1259 1418178901 : return lx==2? -(long)HIGHEXPOBIT: bit_accuracy(lx)-(long)bfffo(*int_MSW(x))-1;
1260 : }
1261 :
1262 : INLINE GEN
1263 140495453 : shiftr(GEN x, long n)
1264 : {
1265 140495453 : const long e = evalexpo(expo(x)+n);
1266 140485947 : const GEN y = rcopy(x);
1267 :
1268 140471886 : if (e & ~EXPOBITS) pari_err_OVERFLOW("expo()");
1269 140473111 : y[1] = (y[1]&~EXPOBITS) | e; return y;
1270 : }
1271 : INLINE GEN
1272 107468879 : mpshift(GEN x,long s) { return (typ(x)==t_INT)?shifti(x,s):shiftr(x,s); }
1273 :
1274 : /* FIXME: adapt/use mpn_[lr]shift instead */
1275 : /* z2[imin..imax] := z1[imin..imax].f shifted left sh bits
1276 : * (feeding f from the right). Assume sh > 0 */
1277 : INLINE void
1278 5897911461 : shift_left(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh)
1279 : {
1280 5897911461 : GEN sb = z1 + imin, se = z1 + imax, te = z2 + imax;
1281 5897911461 : ulong l, m = BITS_IN_LONG - sh, k = f >> m;
1282 40081920057 : while (se > sb) {
1283 34184008596 : l = *se--;
1284 34184008596 : *te-- = (l << sh) | k;
1285 34184008596 : k = l >> m;
1286 : }
1287 5897911461 : *te = (((ulong)*se) << sh) | k;
1288 5897911461 : }
1289 : /* z2[imin..imax] := f.z1[imin..imax-1] shifted right sh bits
1290 : * (feeding f from the left). Assume sh > 0 */
1291 : INLINE void
1292 4439485486 : shift_right(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh)
1293 : {
1294 4439485486 : GEN sb = z1 + imin, se = z1 + imax, tb = z2 + imin;
1295 4439485486 : ulong k, l = *sb++, m = BITS_IN_LONG - sh;
1296 4439485486 : *tb++ = (l >> sh) | (f << m);
1297 23691531993 : while (sb < se) {
1298 19252046507 : k = l << m;
1299 19252046507 : l = *sb++;
1300 19252046507 : *tb++ = (l >> sh) | k;
1301 : }
1302 4439485486 : }
1303 :
1304 : /* Backward compatibility. Inefficient && unused */
1305 : extern ulong hiremainder;
1306 : INLINE ulong
1307 0 : shiftl(ulong x, ulong y)
1308 0 : { hiremainder = x>>(BITS_IN_LONG-y); return (x<<y); }
1309 :
1310 : INLINE ulong
1311 0 : shiftlr(ulong x, ulong y)
1312 0 : { hiremainder = x<<(BITS_IN_LONG-y); return (x>>y); }
1313 :
1314 : INLINE void
1315 302329201 : shiftr_inplace(GEN z, long d)
1316 : {
1317 302329201 : setexpo(z, expo(z)+d);
1318 302327427 : }
1319 :
1320 : /*******************************************************************/
1321 : /* */
1322 : /* ASSIGNMENT */
1323 : /* */
1324 : /*******************************************************************/
1325 : INLINE void
1326 925652996 : affii(GEN x, GEN y)
1327 : {
1328 925652996 : long lx = lgefint(x);
1329 925652996 : if (lg(y)<lx) pari_err_OVERFLOW("t_INT-->t_INT assignment");
1330 37468940571 : while (--lx) y[lx] = x[lx];
1331 925659920 : }
1332 : INLINE void
1333 4049447 : affsi(long s, GEN x)
1334 : {
1335 4049447 : if (!s) x[1] = evalsigne(0) | evallgefint(2);
1336 : else
1337 : {
1338 3906577 : if (s > 0) { x[1] = evalsigne( 1) | evallgefint(3); x[2] = s; }
1339 1134770 : else { x[1] = evalsigne(-1) | evallgefint(3); x[2] = -s; }
1340 : }
1341 4049447 : }
1342 : INLINE void
1343 45103932 : affui(ulong u, GEN x)
1344 : {
1345 45103932 : if (!u) x[1] = evalsigne(0) | evallgefint(2);
1346 45064746 : else { x[1] = evalsigne(1) | evallgefint(3); x[2] = u; }
1347 45103932 : }
1348 :
1349 : INLINE void
1350 419312983 : affsr(long x, GEN y)
1351 : {
1352 419312983 : long sh, i, ly = lg(y);
1353 :
1354 419312983 : if (!x)
1355 : {
1356 70966 : y[1] = evalexpo(-prec2nbits(ly));
1357 70966 : return;
1358 : }
1359 419242017 : if (x < 0) {
1360 3073 : x = -x; sh = bfffo(x);
1361 3073 : y[1] = evalsigne(-1) | _evalexpo((BITS_IN_LONG-1)-sh);
1362 : }
1363 : else
1364 : {
1365 419238944 : sh = bfffo(x);
1366 419238944 : y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
1367 : }
1368 4421094558 : y[2] = ((ulong)x)<<sh; for (i=3; i<ly; i++) y[i]=0;
1369 : }
1370 :
1371 : INLINE void
1372 12412074 : affur(ulong x, GEN y)
1373 : {
1374 12412074 : long sh, i, ly = lg(y);
1375 :
1376 12412074 : if (!x)
1377 : {
1378 1299414 : y[1] = evalexpo(-prec2nbits(ly));
1379 1299414 : return;
1380 : }
1381 11112660 : sh = bfffo(x);
1382 11112660 : y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
1383 49507199 : y[2] = x<<sh; for (i=3; i<ly; i++) y[i] = 0;
1384 : }
1385 :
1386 : INLINE void
1387 266991 : affiz(GEN x, GEN y) { if (typ(y)==t_INT) affii(x,y); else affir(x,y); }
1388 : INLINE void
1389 0 : affsz(long x, GEN y) { if (typ(y)==t_INT) affsi(x,y); else affsr(x,y); }
1390 : INLINE void
1391 658553 : mpaff(GEN x, GEN y) { if (typ(x)==t_INT) affiz(x, y); else affrr(x,y); }
1392 :
1393 : /*******************************************************************/
1394 : /* */
1395 : /* OPERATION + ASSIGNMENT */
1396 : /* */
1397 : /*******************************************************************/
1398 :
1399 0 : INLINE void addiiz(GEN x, GEN y, GEN z)
1400 0 : { pari_sp av = avma; affii(addii(x,y),z); set_avma(av); }
1401 0 : INLINE void addirz(GEN x, GEN y, GEN z)
1402 0 : { pari_sp av = avma; affrr(addir(x,y),z); set_avma(av); }
1403 0 : INLINE void addriz(GEN x, GEN y, GEN z)
1404 0 : { pari_sp av = avma; affrr(addri(x,y),z); set_avma(av); }
1405 1307078 : INLINE void addrrz(GEN x, GEN y, GEN z)
1406 1307078 : { pari_sp av = avma; affrr(addrr(x,y),z); set_avma(av); }
1407 0 : INLINE void addsiz(long s, GEN y, GEN z)
1408 0 : { pari_sp av = avma; affii(addsi(s,y),z); set_avma(av); }
1409 0 : INLINE void addsrz(long s, GEN y, GEN z)
1410 0 : { pari_sp av = avma; affrr(addsr(s,y),z); set_avma(av); }
1411 0 : INLINE void addssz(long s, long y, GEN z)
1412 0 : { pari_sp av = avma; affii(addss(s,y),z); set_avma(av); }
1413 :
1414 0 : INLINE void diviiz(GEN x, GEN y, GEN z)
1415 0 : { pari_sp av = avma; affii(divii(x,y),z); set_avma(av); }
1416 0 : INLINE void divirz(GEN x, GEN y, GEN z)
1417 0 : { pari_sp av = avma; mpaff(divir(x,y),z); set_avma(av); }
1418 0 : INLINE void divisz(GEN x, long y, GEN z)
1419 0 : { pari_sp av = avma; affii(divis(x,y),z); set_avma(av); }
1420 0 : INLINE void divriz(GEN x, GEN y, GEN z)
1421 0 : { pari_sp av = avma; affrr(divri(x,y),z); set_avma(av); }
1422 519 : INLINE void divrrz(GEN x, GEN y, GEN z)
1423 519 : { pari_sp av = avma; affrr(divrr(x,y),z); set_avma(av); }
1424 0 : INLINE void divrsz(GEN y, long s, GEN z)
1425 0 : { pari_sp av = avma; affrr(divrs(y,s),z); set_avma(av); }
1426 0 : INLINE void divsiz(long x, GEN y, GEN z)
1427 0 : { long junk; affsi(sdivsi_rem(x,y,&junk), z); }
1428 0 : INLINE void divsrz(long s, GEN y, GEN z)
1429 0 : { pari_sp av = avma; mpaff(divsr(s,y),z); set_avma(av); }
1430 0 : INLINE void divssz(long x, long y, GEN z)
1431 0 : { affsi(x/y, z); }
1432 :
1433 0 : INLINE void modisz(GEN y, long s, GEN z)
1434 0 : { affsi(smodis(y,s),z); }
1435 0 : INLINE void modsiz(long s, GEN y, GEN z)
1436 0 : { pari_sp av = avma; affii(modsi(s,y),z); set_avma(av); }
1437 0 : INLINE void modssz(long s, long y, GEN z)
1438 0 : { affsi(smodss(s,y),z); }
1439 :
1440 0 : INLINE void mpaddz(GEN x, GEN y, GEN z)
1441 0 : { pari_sp av = avma; mpaff(mpadd(x,y),z); set_avma(av); }
1442 0 : INLINE void mpsubz(GEN x, GEN y, GEN z)
1443 0 : { pari_sp av = avma; mpaff(mpsub(x,y),z); set_avma(av); }
1444 0 : INLINE void mpmulz(GEN x, GEN y, GEN z)
1445 0 : { pari_sp av = avma; mpaff(mpmul(x,y),z); set_avma(av); }
1446 :
1447 0 : INLINE void muliiz(GEN x, GEN y, GEN z)
1448 0 : { pari_sp av = avma; affii(mulii(x,y),z); set_avma(av); }
1449 0 : INLINE void mulirz(GEN x, GEN y, GEN z)
1450 0 : { pari_sp av = avma; mpaff(mulir(x,y),z); set_avma(av); }
1451 0 : INLINE void mulriz(GEN x, GEN y, GEN z)
1452 0 : { pari_sp av = avma; mpaff(mulri(x,y),z); set_avma(av); }
1453 188090 : INLINE void mulrrz(GEN x, GEN y, GEN z)
1454 188090 : { pari_sp av = avma; affrr(mulrr(x,y),z); set_avma(av); }
1455 0 : INLINE void mulsiz(long s, GEN y, GEN z)
1456 0 : { pari_sp av = avma; affii(mulsi(s,y),z); set_avma(av); }
1457 0 : INLINE void mulsrz(long s, GEN y, GEN z)
1458 0 : { pari_sp av = avma; mpaff(mulsr(s,y),z); set_avma(av); }
1459 0 : INLINE void mulssz(long s, long y, GEN z)
1460 0 : { pari_sp av = avma; affii(mulss(s,y),z); set_avma(av); }
1461 :
1462 0 : INLINE void remiiz(GEN x, GEN y, GEN z)
1463 0 : { pari_sp av = avma; affii(remii(x,y),z); set_avma(av); }
1464 0 : INLINE void remisz(GEN y, long s, GEN z)
1465 0 : { pari_sp av = avma; affii(remis(y,s),z); set_avma(av); }
1466 0 : INLINE void remsiz(long s, GEN y, GEN z)
1467 0 : { pari_sp av = avma; affii(remsi(s,y),z); set_avma(av); }
1468 0 : INLINE void remssz(long s, long y, GEN z)
1469 0 : { pari_sp av = avma; affii(remss(s,y),z); set_avma(av); }
1470 :
1471 28 : INLINE void subiiz(GEN x, GEN y, GEN z)
1472 28 : { pari_sp av = avma; affii(subii(x,y),z); set_avma(av); }
1473 0 : INLINE void subirz(GEN x, GEN y, GEN z)
1474 0 : { pari_sp av = avma; affrr(subir(x,y),z); set_avma(av); }
1475 0 : INLINE void subisz(GEN y, long s, GEN z)
1476 0 : { pari_sp av = avma; affii(addsi(-s,y),z); set_avma(av); }
1477 0 : INLINE void subriz(GEN x, GEN y, GEN z)
1478 0 : { pari_sp av = avma; affrr(subri(x,y),z); set_avma(av); }
1479 1296706 : INLINE void subrrz(GEN x, GEN y, GEN z)
1480 1296706 : { pari_sp av = avma; affrr(subrr(x,y),z); set_avma(av); }
1481 0 : INLINE void subrsz(GEN y, long s, GEN z)
1482 0 : { pari_sp av = avma; affrr(addsr(-s,y),z); set_avma(av); }
1483 0 : INLINE void subsiz(long s, GEN y, GEN z)
1484 0 : { pari_sp av = avma; affii(subsi(s,y),z); set_avma(av); }
1485 0 : INLINE void subsrz(long s, GEN y, GEN z)
1486 0 : { pari_sp av = avma; affrr(subsr(s,y),z); set_avma(av); }
1487 0 : INLINE void subssz(long x, long y, GEN z) { addssz(x,-y,z); }
1488 :
1489 : INLINE void
1490 0 : dvmdssz(long x, long y, GEN z, GEN t) {
1491 0 : pari_sp av = avma;
1492 : long r;
1493 0 : affii(divss_rem(x,y, &r), z); set_avma(av); affsi(r,t);
1494 0 : }
1495 : INLINE void
1496 0 : dvmdsiz(long x, GEN y, GEN z, GEN t) {
1497 0 : pari_sp av = avma;
1498 : long r;
1499 0 : affii(divsi_rem(x,y, &r), z); set_avma(av); affsi(r,t);
1500 0 : }
1501 : INLINE void
1502 0 : dvmdisz(GEN x, long y, GEN z, GEN t) {
1503 0 : pari_sp av = avma;
1504 : long r;
1505 0 : affii(divis_rem(x,y, &r),z); set_avma(av); affsi(r,t);
1506 0 : }
1507 : INLINE void
1508 0 : dvmdiiz(GEN x, GEN y, GEN z, GEN t) {
1509 0 : pari_sp av = avma;
1510 : GEN r;
1511 0 : affii(dvmdii(x,y,&r),z); affii(r,t); set_avma(av);
1512 0 : }
|
