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 82745212396 : evallg(long x)
22 : {
23 82745212396 : if (x & ~LGBITS) pari_err_OVERFLOW("lg()");
24 82750157644 : return _evallg(x);
25 : }
26 : INLINE long
27 36273326 : evalvalp(long x)
28 : {
29 36273326 : long v = _evalvalp(x);
30 36273326 : if (v & ~VALPBITS) pari_err_OVERFLOW("valp()");
31 36273396 : return v;
32 : }
33 : INLINE long
34 21316820 : evalvalser(long x)
35 : {
36 21316820 : long v = _evalvalser(x);
37 21316820 : if (v & ~VALSERBITS) pari_err_OVERFLOW("valser()");
38 21316820 : return v;
39 : }
40 : INLINE long
41 10933207496 : evalexpo(long x)
42 : {
43 10933207496 : long v = _evalexpo(x);
44 10933207496 : if (v & ~EXPOBITS) pari_err_OVERFLOW("expo()");
45 10935646614 : return v;
46 : }
47 : INLINE long
48 19354157 : evalprecp(long x)
49 : {
50 19354157 : long v = _evalprecp(x);
51 19354157 : if (x & ~((1UL<<(BITS_IN_LONG-VALPnumBITS))-1)) pari_err_OVERFLOW("precp()");
52 19354178 : return v;
53 : }
54 :
55 : INLINE int
56 161906958 : varncmp(long x, long y)
57 : {
58 161906958 : if (varpriority[x] < varpriority[y]) return 1;
59 125694384 : if (varpriority[x] > varpriority[y]) return -1;
60 68394325 : 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 3132460726 : stackdummy(pari_sp av, pari_sp ltop) {
78 3132460726 : long l = ((GEN)av) - ((GEN)ltop);
79 3132460726 : if (l > 0) {
80 1064812692 : GEN z = (GEN)ltop;
81 1064812692 : 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 3133211399 : }
87 : INLINE void
88 89694310 : fixlg(GEN x, long ly) {
89 89694310 : long lx = lg(x), l = lx - ly;
90 89694310 : if (l > 0)
91 : { /* stackdummy(x+lx, x+ly) */
92 47215451 : GEN z = x + ly;
93 47215451 : z[0] = evaltyp(t_VECSMALL) | evallg(l);
94 47215496 : setlg(x, ly);
95 : #ifdef DEBUG
96 : { long i; for (i = 1; i < l; i++) z[i] = 0; }
97 : #endif
98 : }
99 89694411 : }
100 : /* update lg(z) before affrr(y, z) [ to cater for precision loss ]*/
101 : INLINE void
102 42491406 : 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 >10998*10^7 : set_avma(ulong av) { avma = av; }
113 :
114 : INLINE double
115 173122022 : gc_double(pari_sp av, double d) { set_avma(av); return d; }
116 : INLINE long
117 219713846 : gc_long(pari_sp av, long s) { set_avma(av); return s; }
118 : INLINE ulong
119 29002014 : gc_ulong(pari_sp av, ulong s) { set_avma(av); return s; }
120 : INLINE int
121 47535912 : gc_bool(pari_sp av, int s) { set_avma(av); return s; }
122 : INLINE int
123 2526127 : gc_int(pari_sp av, int s) { set_avma(av); return s; }
124 : INLINE GEN
125 6947022 : gc_NULL(pari_sp av) { set_avma(av); return NULL; }
126 : INLINE GEN
127 13688606184 : gc_const(pari_sp av, GEN x) { set_avma(av); return x; }
128 : INLINE GEN
129 150826 : gc_stoi(pari_sp av, long x) { set_avma(av); return stoi(x); }
130 : INLINE GEN
131 422896 : gc_utoi(pari_sp av, ulong x) { set_avma(av); return utoi(x); }
132 : INLINE GEN
133 1115617 : gc_utoipos(pari_sp av, ulong x) { set_avma(av); return utoipos(x); }
134 :
135 : INLINE GEN
136 79342175549 : new_chunk(size_t x) /* x is a number of longs */
137 : {
138 79342175549 : GEN z = ((GEN) avma) - x;
139 : CHECK_CTRLC
140 79342175549 : if (x > (avma-pari_mainstack->bot) / sizeof(long))
141 14 : new_chunk_resize(x);
142 79342175535 : 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 79310773169 : return z;
155 : }
156 :
157 : INLINE char *
158 10646354 : stack_malloc(size_t N)
159 : {
160 10646354 : long n = nchar2nlong(N);
161 10646341 : return (char*)new_chunk(n);
162 : }
163 :
164 : INLINE char *
165 52195741 : stack_malloc_align(size_t N, long k)
166 : {
167 52195741 : ulong d = ((ulong)avma) % k, e = ((ulong)N) % k;
168 52195741 : if (d) (void)new_chunk(d/sizeof(long));
169 52195742 : if (e) N += k-e;
170 52195742 : return (char*) new_chunk(nchar2nlong(N));
171 : }
172 :
173 : INLINE char *
174 91786 : stack_calloc(size_t N)
175 : {
176 91786 : char *p = stack_malloc(N);
177 91785 : memset(p, 0, N); return p;
178 : }
179 :
180 : INLINE char *
181 956 : stack_calloc_align(size_t N, long k)
182 : {
183 956 : ulong d = ((ulong)avma) % k, e = ((ulong)N) % k;
184 956 : if (d) (void)new_chunk(d/sizeof(long));
185 956 : if (e) N += k-e;
186 956 : return stack_calloc(N);
187 : }
188 :
189 : /* cgetg(lg(x), typ(x)), set *lx. Implicit unsetisclone() */
190 : INLINE GEN
191 1086553598 : cgetg_copy(GEN x, long *plx) {
192 : GEN y;
193 1086553598 : *plx = lg(x); y = new_chunk((size_t)*plx);
194 1086575073 : y[0] = x[0] & (TYPBITS|LGBITS); return y;
195 : }
196 : INLINE GEN
197 393522 : cgetg_block(long x, long y)
198 : {
199 393522 : GEN z = newblock((size_t)x);
200 393410 : z[0] = CLONEBIT | evaltyp(y) | evallg(x);
201 393351 : return z;
202 : }
203 : INLINE GEN
204 21120412051 : cgetg(long x, long y)
205 : {
206 21120412051 : GEN z = new_chunk((size_t)x);
207 21108241941 : z[0] = evaltyp(y) | evallg(x);
208 21093717969 : return z;
209 : }
210 : INLINE GEN
211 23583749873 : cgeti(long x)
212 : {
213 23583749873 : GEN z = new_chunk((size_t)x);
214 23535664705 : z[0] = evaltyp(t_INT) | evallg(x);
215 23517310934 : return z;
216 : }
217 : INLINE GEN
218 13906006204 : cgetipos(long x)
219 : {
220 13906006204 : GEN z = cgeti(x);
221 13882213120 : z[1] = evalsigne(1) | evallgefint(x);
222 13882213120 : return z;
223 : }
224 : INLINE GEN
225 244747198 : cgetineg(long x)
226 : {
227 244747198 : GEN z = cgeti(x);
228 244748941 : z[1] = evalsigne(-1) | evallgefint(x);
229 244748941 : return z;
230 : }
231 : INLINE GEN
232 39653 : cgetr_block(long x)
233 : {
234 39653 : long l = nbits2lg(x);
235 39651 : GEN z = newblock((size_t)l);
236 39663 : z[0] = CLONEBIT | evaltyp(t_REAL) | evallg(l);
237 39659 : return z;
238 : }
239 : INLINE GEN
240 1513213060 : cgetr(long x)
241 : {
242 1513213060 : long l = nbits2lg(x);
243 1513127118 : GEN z = new_chunk((size_t)l);
244 1512616950 : z[0] = evaltyp(t_REAL) | evallg(l);
245 1512306186 : return z;
246 : }
247 :
248 : /*******************************************************************/
249 : /* */
250 : /* COPY, NEGATION, ABSOLUTE VALUE */
251 : /* */
252 : /*******************************************************************/
253 : /* cannot do memcpy because sometimes x and y overlap */
254 : INLINE GEN
255 3333077171 : leafcopy(GEN x)
256 : {
257 3333077171 : long lx = lg(x);
258 3333077171 : GEN y = new_chunk(lx); /* can't use cgetg_copy, in case x,y overlap */
259 17508155831 : while (--lx > 0) y[lx] = x[lx];
260 3332630704 : y[0] = x[0] & (TYPBITS|LGBITS); return y;
261 : }
262 : INLINE GEN
263 7395244509 : icopy(GEN x)
264 : {
265 7395244509 : long i = lgefint(x), lx = i;
266 7395244509 : GEN y = new_chunk(lx); /* can't use cgeti, in case x,y overlap */
267 33303299798 : while (--i > 0) y[i] = x[i];
268 7395186792 : y[0] = evaltyp(t_INT) | evallg(lx);
269 7404351606 : return y;
270 : }
271 : INLINE GEN
272 113969008 : icopyspec(GEN x, long nx)
273 : {
274 113969008 : long i = nx+2, lx = i;
275 113969008 : GEN y = new_chunk(lx); /* can't use cgeti, in case x,y overlap */
276 3057856447 : x -= 2; while (--i >= 2) y[i] = x[i];
277 113968435 : y[1] = evalsigne(1) | evallgefint(lx);
278 113968435 : y[0] = evaltyp(t_INT) | evallg(lx);
279 113968211 : return y;
280 : }
281 844033631 : INLINE GEN rcopy(GEN x) { return leafcopy(x); }
282 707 : INLINE GEN mpcopy(GEN x) { return leafcopy(x); }
283 :
284 : INLINE GEN
285 718896399 : mpabs(GEN x) { GEN y = leafcopy(x); setabssign(y); return y; }
286 : INLINE GEN
287 13410998 : mpabs_shallow(GEN x) { return signe(x) < 0? mpabs(x): x; }
288 659895752 : INLINE GEN absi(GEN x) { return mpabs(x); }
289 40935187 : INLINE GEN absi_shallow(GEN x) { return signe(x) < 0? negi(x): x; }
290 140 : INLINE GEN absr(GEN x) { return mpabs(x); }
291 :
292 : INLINE GEN
293 820498277 : mpneg(GEN x) { GEN y = leafcopy(x); togglesign(y); return y; }
294 573318985 : INLINE GEN negi(GEN x) { return mpneg(x); }
295 2512443 : INLINE GEN negr(GEN x) { return mpneg(x); }
296 :
297 : /* negate in place */
298 : INLINE void
299 1751651203 : togglesign(GEN x) { if (x[1] & SIGNBITS) { x[1] ^= HIGHBIT; } }
300 : INLINE void
301 785995267 : setabssign(GEN x) { x[1] &= ~HIGHBIT; }
302 : /* negate in place, except universal constants */
303 : INLINE void
304 122190654 : togglesign_safe(GEN *px)
305 : {
306 122190654 : switch(*px - gen_1) /* gen_1, gen_2, gen_m1, gen_m2 */
307 : {
308 2472518 : case 0: *px = gen_m1; break;
309 4 : case 3: *px = gen_m2; break;
310 557211 : case 6: *px = gen_1; break;
311 0 : case 9: *px = gen_2; break;
312 119160921 : default: togglesign(*px);
313 : }
314 122192882 : }
315 : /* setsigne(y, signe(x)) */
316 : INLINE void
317 0 : affectsign(GEN x, GEN y)
318 : {
319 0 : y[1] = (x[1] & SIGNBITS) | (y[1] & ~SIGNBITS);
320 0 : }
321 : /* copies sign in place, except for universal constants */
322 : INLINE void
323 9994817 : affectsign_safe(GEN x, GEN *py)
324 : {
325 9994817 : if (((*py)[1] ^ x[1]) & HIGHBIT) togglesign_safe(py);
326 9994817 : }
327 : /*******************************************************************/
328 : /* */
329 : /* GEN -> LONG, LONG -> GEN */
330 : /* */
331 : /*******************************************************************/
332 : /* assume x != 0, return -x as a t_INT */
333 : INLINE GEN
334 243898760 : utoineg(ulong x) { GEN y = cgetineg(3); y[2] = x; return y; }
335 : /* assume x != 0, return utoi(x) */
336 : INLINE GEN
337 12082194951 : utoipos(ulong x) { GEN y = cgetipos(3); y[2] = x; return y; }
338 : INLINE GEN
339 10215628128 : utoi(ulong x) { return x? utoipos(x): gen_0; }
340 : INLINE GEN
341 713358566 : stoi(long x)
342 : {
343 713358566 : if (!x) return gen_0;
344 494650908 : return x > 0? utoipos((ulong)x): utoineg((ulong)-x);
345 : }
346 :
347 : /* x 2^BIL + y */
348 : INLINE GEN
349 7541996754 : uutoi(ulong x, ulong y)
350 : {
351 : GEN z;
352 7541996754 : if (!x) return utoi(y);
353 724598802 : z = cgetipos(4);
354 727089748 : *int_W_lg(z, 1, 4) = x;
355 727089748 : *int_W_lg(z, 0, 4) = y; return z;
356 : }
357 : /* - (x 2^BIL + y) */
358 : INLINE GEN
359 246156 : uutoineg(ulong x, ulong y)
360 : {
361 : GEN z;
362 246156 : if (!x) return y? utoineg(y): gen_0;
363 10425 : z = cgetineg(4);
364 10425 : *int_W_lg(z, 1, 4) = x;
365 10425 : *int_W_lg(z, 0, 4) = y; return z;
366 : }
367 :
368 : INLINE long
369 427419022 : itos(GEN x)
370 : {
371 427419022 : long s = signe(x);
372 : long u;
373 :
374 427419022 : if (!s) return 0;
375 400440388 : u = x[2];
376 400440388 : if (lgefint(x) > 3 || u < 0)
377 27 : pari_err_OVERFLOW("t_INT-->long assignment");
378 400438317 : return (s>0) ? u : -u;
379 : }
380 : /* as itos, but return 0 if too large. Cf is_bigint */
381 : INLINE long
382 17526756 : itos_or_0(GEN x) {
383 : long n;
384 17526756 : if (lgefint(x) != 3 || (n = x[2]) & HIGHBIT) return 0;
385 15782753 : return signe(x) > 0? n: -n;
386 : }
387 : INLINE ulong
388 151992306 : itou(GEN x)
389 : {
390 151992306 : switch(lgefint(x)) {
391 10044229 : case 2: return 0;
392 141948178 : case 3: return x[2];
393 0 : default:
394 0 : pari_err_OVERFLOW("t_INT-->ulong assignment");
395 : return 0; /* LCOV_EXCL_LINE */
396 : }
397 : }
398 :
399 : /* as itou, but return 0 if too large. Cf is_bigint */
400 : INLINE ulong
401 2469097 : itou_or_0(GEN x) {
402 2469097 : if (lgefint(x) != 3) return 0;
403 2448026 : return (ulong)x[2];
404 : }
405 :
406 : INLINE ulong
407 5520653 : umuluu_or_0(ulong x, ulong y)
408 : {
409 : ulong z;
410 : LOCAL_HIREMAINDER;
411 5520653 : z = mulll(x, y);
412 5520653 : return hiremainder? 0: z;
413 : }
414 : /* return x*y if <= n, else 0. Beware overflow */
415 : INLINE ulong
416 5613979 : umuluu_le(ulong x, ulong y, ulong n)
417 : {
418 : ulong z;
419 : LOCAL_HIREMAINDER;
420 5613979 : z = mulll(x, y);
421 5613979 : return (hiremainder || z > n)? 0: z;
422 : }
423 :
424 : INLINE GEN
425 184218694 : real_0_bit(long bitprec) { GEN x=cgetg(2, t_REAL); x[1]=evalexpo(bitprec); return x; }
426 : INLINE GEN
427 629488 : real_0(long prec) { return real_0_bit(-prec2nbits(prec)); }
428 : INLINE GEN
429 2675822 : real_1_bit(long bit) { return real_1(nbits2prec(bit)); }
430 : INLINE GEN
431 104795035 : real_1(long prec) {
432 104795035 : GEN x = cgetr(prec);
433 104786531 : long i, l = lg(x);
434 104786531 : x[1] = evalsigne(1) | _evalexpo(0);
435 474667719 : x[2] = (long)HIGHBIT; for (i=3; i<l; i++) x[i] = 0;
436 104786531 : return x;
437 : }
438 : INLINE GEN
439 329 : real_m1(long prec) {
440 329 : GEN x = cgetr(prec);
441 329 : long i, l = lg(x);
442 329 : x[1] = evalsigne(-1) | _evalexpo(0);
443 1617 : x[2] = (long)HIGHBIT; for (i=3; i<l; i++) x[i] = 0;
444 329 : return x;
445 : }
446 :
447 : /* 2.^n */
448 : INLINE GEN
449 591162 : real2n(long n, long prec) { GEN z = real_1(prec); setexpo(z, n); return z; }
450 : INLINE GEN
451 0 : real_m2n(long n, long prec) { GEN z = real_m1(prec); setexpo(z, n); return z; }
452 : INLINE GEN
453 413958156 : stor(long s, long prec) { GEN z = cgetr(prec); affsr(s,z); return z; }
454 : INLINE GEN
455 12213273 : utor(ulong s, long prec){ GEN z = cgetr(prec); affur(s,z); return z; }
456 : INLINE GEN
457 650697319 : itor(GEN x, long prec) { GEN z = cgetr(prec); affir(x,z); return z; }
458 : INLINE GEN
459 230208381 : rtor(GEN x, long prec) { GEN z = cgetr(prec); affrr(x,z); return z; }
460 :
461 22121536 : INLINE ulong int_bit(GEN x, long n)
462 : {
463 22121536 : long r, q = dvmdsBIL(n, &r);
464 22109525 : return q < lgefint(x)-2?((ulong)*int_W(x,q) >> r) & 1UL:0;
465 : }
466 :
467 : /*******************************************************************/
468 : /* */
469 : /* COMPARISON */
470 : /* */
471 : /*******************************************************************/
472 : INLINE int
473 1289805 : cmpss(long a, long b)
474 1289805 : { return a>b? 1: (a<b? -1: 0); }
475 :
476 : INLINE int
477 1418182584 : cmpuu(ulong a, ulong b)
478 1418182584 : { return a>b? 1: (a<b? -1: 0); }
479 :
480 : INLINE int
481 6190797 : cmpir(GEN x, GEN y)
482 : {
483 : pari_sp av;
484 : GEN z;
485 :
486 6190797 : if (!signe(x)) return -signe(y);
487 477245 : if (!signe(y))
488 : {
489 2248 : if (expo(y) >= expi(x)) return 0;
490 2220 : return signe(x);
491 : }
492 474997 : av=avma; z = itor(x, realprec(y)); set_avma(av);
493 474987 : return cmprr(z,y); /* cmprr does no memory adjustment */
494 : }
495 : INLINE int
496 408913 : cmpri(GEN x, GEN y) { return -cmpir(y,x); }
497 : INLINE int
498 586514 : cmpsr(long x, GEN y)
499 : {
500 : pari_sp av;
501 : GEN z;
502 :
503 586514 : if (!x) return -signe(y);
504 586514 : av=avma; z = stor(x, LOWDEFAULTPREC); set_avma(av);
505 586514 : return cmprr(z,y);
506 : }
507 : INLINE int
508 40996 : cmprs(GEN x, long y) { return -cmpsr(y,x); }
509 : /* compare x and y */
510 : INLINE int
511 8720969 : cmpui(ulong x, GEN y)
512 : {
513 : ulong p;
514 8720969 : if (!x) return -signe(y);
515 8720906 : if (signe(y) <= 0) return 1;
516 8609137 : if (lgefint(y) > 3) return -1;
517 8380939 : p = y[2]; if (p == x) return 0;
518 8310942 : return p < x ? 1 : -1;
519 : }
520 : INLINE int
521 8716732 : cmpiu(GEN x, ulong y) { return -cmpui(y,x); }
522 : /* compare x and |y| */
523 : INLINE int
524 32013873 : abscmpui(ulong x, GEN y)
525 : {
526 32013873 : long l = lgefint(y);
527 : ulong p;
528 :
529 32013873 : if (!x) return (l > 2)? -1: 0;
530 32013859 : if (l == 2) return 1;
531 31666730 : if (l > 3) return -1;
532 31638667 : p = y[2]; if (p == x) return 0;
533 31057482 : return p < x ? 1 : -1;
534 : }
535 : INLINE int
536 32012398 : abscmpiu(GEN x, ulong y) { return -abscmpui(y,x); }
537 : INLINE int
538 3377446 : cmpsi(long x, GEN y)
539 : {
540 : ulong p;
541 :
542 3377446 : if (!x) return -signe(y);
543 :
544 3374274 : if (x > 0)
545 : {
546 3373294 : if (signe(y)<=0) return 1;
547 3373000 : if (lgefint(y)>3) return -1;
548 3357692 : p = y[2]; if (p == (ulong)x) return 0;
549 3287962 : return p < (ulong)x ? 1 : -1;
550 : }
551 :
552 980 : if (signe(y)>=0) return -1;
553 119 : if (lgefint(y)>3) return 1;
554 119 : p = y[2]; if (p == (ulong)-x) return 0;
555 14 : return p < (ulong)(-x) ? -1 : 1;
556 : }
557 : INLINE int
558 3352335 : cmpis(GEN x, long y) { return -cmpsi(y,x); }
559 : INLINE int
560 2126083 : mpcmp(GEN x, GEN y)
561 : {
562 2126083 : if (typ(x)==t_INT)
563 60304 : return (typ(y)==t_INT) ? cmpii(x,y) : cmpir(x,y);
564 2065779 : return (typ(y)==t_INT) ? -cmpir(y,x) : cmprr(x,y);
565 : }
566 :
567 : /* x == y ? */
568 : INLINE int
569 2917094 : equalui(ulong x, GEN y)
570 : {
571 2917094 : if (!x) return !signe(y);
572 2916401 : if (signe(y) <= 0 || lgefint(y) != 3) return 0;
573 2904837 : return ((ulong)y[2] == (ulong)x);
574 : }
575 : /* x == y ? */
576 : INLINE int
577 524098 : equalsi(long x, GEN y)
578 : {
579 524098 : if (!x) return !signe(y);
580 524098 : if (x > 0)
581 : {
582 521899 : if (signe(y) <= 0 || lgefint(y) != 3) return 0;
583 469449 : return ((ulong)y[2] == (ulong)x);
584 : }
585 2199 : if (signe(y) >= 0 || lgefint(y) != 3) return 0;
586 2078 : return ((ulong)y[2] == (ulong)-x);
587 : }
588 : /* x == |y| ? */
589 : INLINE int
590 38006274 : absequalui(ulong x, GEN y)
591 : {
592 38006274 : if (!x) return !signe(y);
593 38006274 : return (lgefint(y) == 3 && (ulong)y[2] == x);
594 : }
595 : INLINE int
596 36311666 : absequaliu(GEN x, ulong y) { return absequalui(y,x); }
597 : INLINE int
598 523916 : equalis(GEN x, long y) { return equalsi(y,x); }
599 : INLINE int
600 2917095 : equaliu(GEN x, ulong y) { return equalui(y,x); }
601 :
602 : /* assume x != 0, is |x| == 2^n ? */
603 : INLINE int
604 1234799 : absrnz_equal2n(GEN x) {
605 1234799 : if ((ulong)x[2]==HIGHBIT)
606 : {
607 29641 : long i, lx = lg(x);
608 92880 : for (i = 3; i < lx; i++)
609 67301 : if (x[i]) return 0;
610 25579 : return 1;
611 : }
612 1205158 : return 0;
613 : }
614 : /* assume x != 0, is |x| == 1 ? */
615 : INLINE int
616 3684791 : absrnz_equal1(GEN x) { return !expo(x) && absrnz_equal2n(x); }
617 :
618 : INLINE long
619 8950233178 : maxss(long x, long y) { return x>y?x:y; }
620 : INLINE long
621 1658360620 : minss(long x, long y) { return x<y?x:y; }
622 : INLINE long
623 7566352 : minuu(ulong x, ulong y) { return x<y?x:y; }
624 : INLINE long
625 14280174 : maxuu(ulong x, ulong y) { return x>y?x:y; }
626 : INLINE double
627 2906588 : maxdd(double x, double y) { return x>y?x:y; }
628 : INLINE double
629 249264 : mindd(double x, double y) { return x<y?x:y; }
630 :
631 : /*******************************************************************/
632 : /* */
633 : /* ADD / SUB */
634 : /* */
635 : /*******************************************************************/
636 : INLINE GEN
637 25067 : subuu(ulong x, ulong y)
638 : {
639 : ulong z;
640 : LOCAL_OVERFLOW;
641 25067 : z = subll(x, y);
642 25067 : return overflow? utoineg(-z): utoi(z);
643 : }
644 : INLINE GEN
645 3064210834 : adduu(ulong x, ulong y) { ulong t = x+y; return uutoi((t < x), t); }
646 :
647 : INLINE GEN
648 25067 : addss(long x, long y)
649 : {
650 25067 : if (!x) return stoi(y);
651 25067 : if (!y) return stoi(x);
652 25067 : if (x > 0) return y > 0? adduu(x,y): subuu(x, -y);
653 :
654 25067 : if (y > 0) return subuu(y, -x);
655 : else { /* - adduu(-x, -y) */
656 0 : ulong t = (-x)+(-y); return uutoineg((t < (ulong)(-x)), t);
657 : }
658 : }
659 25067 : INLINE GEN subss(long x, long y) { return addss(-y,x); }
660 :
661 : INLINE GEN
662 7413571976 : subii(GEN x, GEN y)
663 : {
664 7413571976 : if (x==y) return gen_0; /* frequent with x = y = gen_0 */
665 6089994981 : return addii_sign(x, signe(x), y, -signe(y));
666 : }
667 : INLINE GEN
668 10111530320 : addii(GEN x, GEN y) { return addii_sign(x, signe(x), y, signe(y)); }
669 : INLINE GEN
670 2426634914 : addrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, signe(y)); }
671 : INLINE GEN
672 873955177 : subrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, -signe(y)); }
673 : INLINE GEN
674 461619936 : addir(GEN x, GEN y) { return addir_sign(x, signe(x), y, signe(y)); }
675 : INLINE GEN
676 2801174 : subir(GEN x, GEN y) { return addir_sign(x, signe(x), y, -signe(y)); }
677 : INLINE GEN
678 5957699 : subri(GEN x, GEN y) { return addir_sign(y, -signe(y), x, signe(x)); }
679 : INLINE GEN
680 270932395 : addsi(long x, GEN y) { return addsi_sign(x, y, signe(y)); }
681 : INLINE GEN
682 91964824 : addui(ulong x, GEN y) { return addui_sign(x, y, signe(y)); }
683 : INLINE GEN
684 5819091 : subsi(long x, GEN y) { return addsi_sign(x, y, -signe(y)); }
685 : INLINE GEN
686 126536399 : subui(ulong x, GEN y) { return addui_sign(x, y, -signe(y)); }
687 :
688 : /*******************************************************************/
689 : /* */
690 : /* MOD, REM, DIV */
691 : /* */
692 : /*******************************************************************/
693 92282551 : INLINE ulong mod2BIL(GEN x) { return *int_LSW(x); }
694 0 : INLINE long mod64(GEN x) { return mod2BIL(x) & 63; }
695 259 : INLINE long mod32(GEN x) { return mod2BIL(x) & 31; }
696 236257 : INLINE long mod16(GEN x) { return mod2BIL(x) & 15; }
697 12779052 : INLINE long mod8(GEN x) { return mod2BIL(x) & 7; }
698 4083931 : INLINE long mod4(GEN x) { return mod2BIL(x) & 3; }
699 52828636 : INLINE long mod2(GEN x) { return mod2BIL(x) & 1; }
700 : INLINE int
701 83887787 : mpodd(GEN x) { return signe(x) && mod2(x); }
702 : /* x mod 2^n, n < BITS_IN_LONG */
703 : INLINE ulong
704 38538490 : umodi2n(GEN x, long n)
705 : {
706 38538490 : long s = signe(x);
707 38538490 : const ulong _2n = 1UL << n;
708 : ulong m;
709 38538490 : if (!s) return 0;
710 37077189 : m = *int_LSW(x) & (_2n - 1);
711 37077189 : if (s < 0 && m) m = _2n - m;
712 37077189 : return m;
713 : }
714 0 : INLINE ulong Mod64(GEN x){ return umodi2n(x,6); }
715 167650 : INLINE ulong Mod32(GEN x){ return umodi2n(x,5); }
716 244826 : INLINE ulong Mod16(GEN x){ return umodi2n(x,4); }
717 1780948 : INLINE ulong Mod8(GEN x) { return umodi2n(x,3); }
718 34655164 : INLINE ulong Mod4(GEN x) { return umodi2n(x,2); }
719 1689254 : INLINE ulong Mod2(GEN x) { return umodi2n(x,1); }
720 :
721 : INLINE GEN
722 45831509 : truedivii(GEN a,GEN b) { return truedvmdii(a,b,NULL); }
723 : INLINE GEN
724 248133 : truedivis(GEN a, long b) { return truedvmdis(a,b,NULL); }
725 : INLINE GEN
726 6197775 : truedivsi(long a, GEN b) { return truedvmdsi(a,b,NULL); }
727 :
728 : INLINE GEN
729 12113747 : divii(GEN a, GEN b) { return dvmdii(a,b,NULL); }
730 : INLINE GEN
731 2454451550 : remii(GEN a, GEN b) { return dvmdii(a,b,ONLY_REM); }
732 :
733 : INLINE GEN
734 0 : divss(long x, long y) { return stoi(x / y); }
735 : INLINE GEN
736 0 : modss(long x, long y) { return utoi(smodss(x, y)); }
737 : INLINE GEN
738 0 : remss(long x, long y) { return stoi(x % y); }
739 : INLINE long
740 6414 : smodss(long x, long y)
741 : {
742 6414 : long r = x%y;
743 6414 : return (r >= 0)? r: labs(y) + r;
744 : }
745 : INLINE ulong
746 716025521 : umodsu(long x, ulong y)
747 : {
748 716025521 : return x>=0 ? x%y: Fl_neg((-x)%y, y);
749 : }
750 :
751 : INLINE long
752 0 : sdivss_rem(long x, long y, long *r)
753 : {
754 : long q;
755 : LOCAL_HIREMAINDER;
756 0 : if (!y) pari_err_INV("sdivss_rem",gen_0);
757 0 : hiremainder = 0; q = divll((ulong)labs(x),(ulong)labs(y));
758 0 : if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
759 0 : if (y < 0) q = -q;
760 0 : *r = hiremainder; return q;
761 : }
762 : INLINE GEN
763 0 : divss_rem(long x, long y, long *r) { return stoi(sdivss_rem(x,y,r)); }
764 : INLINE ulong
765 158497970 : udivuu_rem(ulong x, ulong y, ulong *r)
766 : {
767 158497970 : if (!y) pari_err_INV("udivuu_rem",gen_0);
768 158497970 : *r = x % y; return x / y;
769 : }
770 : INLINE ulong
771 2012318 : ceildivuu(ulong a, ulong b)
772 : {
773 : ulong c;
774 2012318 : if (!a) return 0;
775 1834039 : c = a / b; return (a % b)? c+1: c;
776 : }
777 :
778 : INLINE ulong
779 13282 : uabsdivui_rem(ulong x, GEN y, ulong *r)
780 : {
781 13282 : long q, s = signe(y);
782 : LOCAL_HIREMAINDER;
783 :
784 13282 : if (!s) pari_err_INV("uabsdivui_rem",gen_0);
785 13282 : if (!x || lgefint(y)>3) { *r = x; return 0; }
786 12618 : hiremainder=0; q = (long)divll(x, (ulong)y[2]);
787 12618 : if (s < 0) q = -q;
788 12618 : *r = hiremainder; return q;
789 : }
790 :
791 : /* assume d != 0 and |n| / d can be represented as an ulong.
792 : * Return |n|/d, set *r = |n| % d */
793 : INLINE ulong
794 8349446 : uabsdiviu_rem(GEN n, ulong d, ulong *r)
795 : {
796 8349446 : switch(lgefint(n))
797 : {
798 0 : case 2: *r = 0; return 0;
799 8349446 : case 3:
800 : {
801 8349446 : ulong nn = n[2];
802 8349446 : *r = nn % d; return nn / d;
803 : }
804 0 : default: /* 4 */
805 : {
806 : ulong n1, n0, q;
807 : LOCAL_HIREMAINDER;
808 0 : n0 = *int_W(n,0);
809 0 : n1 = *int_W(n,1);
810 0 : hiremainder = n1;
811 0 : q = divll(n0, d);
812 0 : *r = hiremainder; return q;
813 : }
814 : }
815 : }
816 :
817 : INLINE long
818 51288087 : sdivsi_rem(long x, GEN y, long *r)
819 : {
820 51288087 : long q, s = signe(y);
821 : LOCAL_HIREMAINDER;
822 :
823 51288087 : if (!s) pari_err_INV("sdivsi_rem",gen_0);
824 51288087 : if (!x || lgefint(y)>3 || ((long)y[2]) < 0) { *r = x; return 0; }
825 49333834 : hiremainder=0; q = (long)divll(labs(x), (ulong)y[2]);
826 49333834 : if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
827 49333834 : if (s < 0) q = -q;
828 49333834 : *r = hiremainder; return q;
829 : }
830 : INLINE GEN
831 0 : divsi_rem(long s, GEN y, long *r) { return stoi(sdivsi_rem(s,y,r)); }
832 :
833 : INLINE long
834 100943 : sdivsi(long x, GEN y)
835 : {
836 100943 : long q, s = signe(y);
837 :
838 100943 : if (!s) pari_err_INV("sdivsi",gen_0);
839 100943 : if (!x || lgefint(y)>3 || ((long)y[2]) < 0) return 0;
840 100838 : q = labs(x) / y[2];
841 100838 : if (x < 0) q = -q;
842 100838 : if (s < 0) q = -q;
843 100838 : return q;
844 : }
845 :
846 : INLINE GEN
847 0 : dvmdss(long x, long y, GEN *z)
848 : {
849 : long r;
850 0 : GEN q = divss_rem(x,y, &r);
851 0 : *z = stoi(r); return q;
852 : }
853 : INLINE long
854 5961345732 : dvmdsBIL(long n, long *r) { *r = remsBIL(n); return divsBIL(n); }
855 : INLINE ulong
856 164569083 : dvmduBIL(ulong n, ulong *r) { *r = remsBIL(n); return divsBIL(n); }
857 : INLINE GEN
858 0 : dvmdsi(long x, GEN y, GEN *z)
859 : {
860 : long r;
861 0 : GEN q = divsi_rem(x,y, &r);
862 0 : *z = stoi(r); return q;
863 : }
864 : INLINE GEN
865 0 : dvmdis(GEN x, long y, GEN *z)
866 : {
867 : long r;
868 0 : GEN q = divis_rem(x,y, &r);
869 0 : *z = stoi(r); return q;
870 : }
871 :
872 : INLINE long
873 21006743 : smodis(GEN x, long y)
874 : {
875 21006743 : pari_sp av = avma;
876 21006743 : long r; (void)divis_rem(x,y, &r);
877 21006743 : return gc_long(av, (r >= 0)? r: labs(y) + r);
878 : }
879 : INLINE GEN
880 19602368 : modis(GEN x, long y) { return stoi(smodis(x,y)); }
881 : INLINE GEN
882 45090067 : modsi(long x, GEN y) {
883 45090067 : long r; (void)sdivsi_rem(x, y, &r);
884 45090067 : return (r >= 0)? stoi(r): addsi_sign(r, y, 1);
885 : }
886 :
887 : INLINE ulong
888 1511727 : umodui(ulong x, GEN y)
889 : {
890 1511727 : if (!signe(y)) pari_err_INV("umodui",gen_0);
891 1511727 : if (!x || lgefint(y) > 3) return x;
892 415241 : return x % (ulong)y[2];
893 : }
894 :
895 : INLINE ulong
896 210587 : ugcdiu(GEN x, ulong y) { return ugcd(umodiu(x,y), y); }
897 : INLINE ulong
898 2737 : ugcdui(ulong y, GEN x) { return ugcd(umodiu(x,y), y); }
899 :
900 : INLINE GEN
901 0 : remsi(long x, GEN y)
902 0 : { long r; (void)sdivsi_rem(x,y, &r); return stoi(r); }
903 : INLINE GEN
904 0 : remis(GEN x, long y)
905 : {
906 0 : pari_sp av = avma;
907 : long r;
908 0 : (void)divis_rem(x,y, &r); return gc_stoi(av, r);
909 : }
910 :
911 : INLINE GEN
912 0 : rdivis(GEN x, long y, long prec)
913 : {
914 0 : GEN z = cgetr(prec);
915 0 : pari_sp av = avma;
916 0 : affrr(divrs(itor(x,prec), y),z);
917 0 : set_avma(av); return z;
918 : }
919 : INLINE GEN
920 0 : rdivsi(long x, GEN y, long prec)
921 : {
922 0 : GEN z = cgetr(prec);
923 0 : pari_sp av = avma;
924 0 : affrr(divsr(x, itor(y,prec)), z);
925 0 : set_avma(av); return z;
926 : }
927 : INLINE GEN
928 839647 : rdivss(long x, long y, long prec)
929 : {
930 839647 : GEN z = cgetr(prec);
931 839647 : pari_sp av = avma;
932 839647 : affrr(divrs(stor(x, prec), y), z);
933 839647 : set_avma(av); return z;
934 : }
935 :
936 : INLINE void
937 7801768 : rdiviiz(GEN x, GEN y, GEN z)
938 : {
939 7801768 : long lz = lg(z), lx = lgefint(x), ly = lgefint(y);
940 7801768 : if (lx == 2) { affur(0, z); return; }
941 7801768 : if (ly == 3)
942 : {
943 2212764 : affir(x, z); if (signe(y) < 0) togglesign(z);
944 2212734 : affrr(divru(z, y[2]), z);
945 : }
946 5589004 : else if (lx > lz + 1 || ly > lz + 1)
947 : {
948 1851873 : affir(x,z); affrr(divri(z, y), z);
949 : }
950 : else
951 : {
952 3737131 : long b = prec2nbits(lg2prec(lz)) + expi(y) - expi(x) + 1;
953 3737191 : GEN q = divii(b > 0? shifti(x, b): x, y);
954 3737198 : affir(q, z); if (b > 0) shiftr_inplace(z, -b);
955 : }
956 7801889 : set_avma((ulong)z);
957 : }
958 : INLINE GEN
959 7762239 : rdivii(GEN x, GEN y, long prec)
960 7762239 : { GEN z = cgetr(prec); rdiviiz(x, y, z); return z; }
961 : INLINE GEN
962 7346555 : fractor(GEN x, long prec)
963 7346555 : { return rdivii(gel(x,1), gel(x,2), prec); }
964 :
965 : INLINE int
966 16089264 : dvdii(GEN x, GEN y)
967 : {
968 16089264 : pari_sp av = avma;
969 : GEN r;
970 16089264 : if (!signe(x)) return 1;
971 14795784 : if (!signe(y)) return 0;
972 14795784 : r = remii(x,y);
973 14804888 : return gc_bool(av, r == gen_0);
974 : }
975 : INLINE int
976 371 : dvdsi(long x, GEN y)
977 : {
978 371 : if (x == 0) return 1;
979 266 : if (!signe(y)) return 0;
980 266 : if (lgefint(y) != 3) return 0;
981 259 : return x % y[2] == 0;
982 : }
983 : INLINE int
984 167195 : dvdui(ulong x, GEN y)
985 : {
986 167195 : if (x == 0) return 1;
987 167195 : if (!signe(y)) return 0;
988 167195 : if (lgefint(y) != 3) return 0;
989 156574 : return x % y[2] == 0;
990 : }
991 : INLINE int
992 33492 : dvdis(GEN x, long y)
993 33492 : { return y? smodis(x, y) == 0: signe(x) == 0; }
994 : INLINE int
995 576505 : dvdiu(GEN x, ulong y)
996 576505 : { return y? umodiu(x, y) == 0: signe(x) == 0; }
997 :
998 : INLINE int
999 0 : dvdisz(GEN x, long y, GEN z)
1000 : {
1001 0 : const pari_sp av = avma;
1002 : long r;
1003 0 : GEN p1 = divis_rem(x,y, &r);
1004 0 : set_avma(av); if (r) return 0;
1005 0 : affii(p1,z); return 1;
1006 : }
1007 : INLINE int
1008 0 : dvdiuz(GEN x, ulong y, GEN z)
1009 : {
1010 0 : const pari_sp av = avma;
1011 : ulong r;
1012 0 : GEN p1 = absdiviu_rem(x,y, &r);
1013 0 : set_avma(av); if (r) return 0;
1014 0 : affii(p1,z); return 1;
1015 : }
1016 : INLINE int
1017 5431 : dvdiiz(GEN x, GEN y, GEN z)
1018 : {
1019 5431 : const pari_sp av=avma;
1020 5431 : GEN p2, p1 = dvmdii(x,y,&p2);
1021 5431 : if (signe(p2)) return gc_bool(av,0);
1022 4230 : affii(p1,z); return gc_bool(av,1);
1023 : }
1024 :
1025 : INLINE ulong
1026 81821485 : remlll_pre(ulong u2, ulong u1, ulong u0, ulong n, ulong ninv)
1027 : {
1028 81821485 : u1 = remll_pre(u2, u1, n, ninv);
1029 82807139 : return remll_pre(u1, u0, n, ninv);
1030 : }
1031 :
1032 : INLINE ulong
1033 1985012256 : Fl_sqr_pre(ulong a, ulong p, ulong pi)
1034 : {
1035 : ulong x;
1036 : LOCAL_HIREMAINDER;
1037 1985012256 : x = mulll(a,a);
1038 1985012256 : return remll_pre(hiremainder, x, p, pi);
1039 : }
1040 :
1041 : INLINE ulong
1042 3318584794 : Fl_mul_pre(ulong a, ulong b, ulong p, ulong pi)
1043 : {
1044 : ulong x;
1045 : LOCAL_HIREMAINDER;
1046 3318584794 : x = mulll(a,b);
1047 3318584794 : return remll_pre(hiremainder, x, p, pi);
1048 : }
1049 :
1050 : INLINE ulong
1051 7051050353 : Fl_addmul_pre(ulong y0, ulong x0, ulong x1, ulong p, ulong pi)
1052 : {
1053 : ulong l0, h0;
1054 : LOCAL_HIREMAINDER;
1055 7051050353 : hiremainder = y0;
1056 7051050353 : l0 = addmul(x0, x1); h0 = hiremainder;
1057 7051050353 : return remll_pre(h0, l0, p, pi);
1058 : }
1059 :
1060 : INLINE ulong
1061 55556211 : Fl_addmulmul_pre(ulong x0, ulong y0, ulong x1, ulong y1, ulong p, ulong pi)
1062 : {
1063 : ulong l0, l1, h0, h1;
1064 : LOCAL_OVERFLOW;
1065 : LOCAL_HIREMAINDER;
1066 55556211 : l0 = mulll(x0, y0); h0 = hiremainder;
1067 55556211 : l1 = mulll(x1, y1); h1 = hiremainder;
1068 55556211 : l0 = addll(l0, l1); h0 = addllx(h0, h1);
1069 55556211 : return overflow ? remlll_pre(1, h0, l0, p, pi): remll_pre(h0, l0, p, pi);
1070 : }
1071 :
1072 : INLINE ulong
1073 219824 : Fl_ellj_pre(ulong a4, ulong a6, ulong p, ulong pi)
1074 : {
1075 : /* a43 = 4 a4^3 */
1076 219824 : ulong a43 = Fl_double(Fl_double(
1077 : Fl_mul_pre(a4, Fl_sqr_pre(a4, p, pi), p, pi), p), p);
1078 : /* a62 = 27 a6^2 */
1079 219829 : ulong a62 = Fl_mul_pre(Fl_sqr_pre(a6, p, pi), 27 % p, p, pi);
1080 219827 : ulong z1 = Fl_mul_pre(a43, 1728 % p, p, pi);
1081 219826 : ulong z2 = Fl_add(a43, a62, p);
1082 219824 : return Fl_div(z1, z2, p);
1083 : }
1084 :
1085 : /*******************************************************************/
1086 : /* */
1087 : /* MP (INT OR REAL) */
1088 : /* */
1089 : /*******************************************************************/
1090 : INLINE GEN
1091 49 : mptrunc(GEN x) { return typ(x)==t_INT? icopy(x): truncr(x); }
1092 : INLINE GEN
1093 0 : mpfloor(GEN x) { return typ(x)==t_INT? icopy(x): floorr(x); }
1094 : INLINE GEN
1095 0 : mpceil(GEN x) { return typ(x)==t_INT? icopy(x): ceilr(x); }
1096 : INLINE GEN
1097 1848110 : mpround(GEN x) { return typ(x) == t_INT? icopy(x): roundr(x); }
1098 :
1099 : INLINE long
1100 32488884 : mpexpo(GEN x) { return typ(x) == t_INT? expi(x): expo(x); }
1101 :
1102 : INLINE GEN
1103 353052048 : mpadd(GEN x, GEN y)
1104 : {
1105 353052048 : if (typ(x)==t_INT)
1106 20173776 : return (typ(y)==t_INT) ? addii(x,y) : addir(x,y);
1107 332878272 : return (typ(y)==t_INT) ? addir(y,x) : addrr(x,y);
1108 : }
1109 : INLINE GEN
1110 229821964 : mpsub(GEN x, GEN y)
1111 : {
1112 229821964 : if (typ(x)==t_INT)
1113 500340 : return (typ(y)==t_INT) ? subii(x,y) : subir(x,y);
1114 229321624 : return (typ(y)==t_INT) ? subri(x,y) : subrr(x,y);
1115 : }
1116 : INLINE GEN
1117 582049751 : mpmul(GEN x, GEN y)
1118 : {
1119 582049751 : if (typ(x)==t_INT)
1120 36644415 : return (typ(y)==t_INT) ? mulii(x,y) : mulir(x,y);
1121 545405336 : return (typ(y)==t_INT) ? mulir(y,x) : mulrr(x,y);
1122 : }
1123 : INLINE GEN
1124 74508277 : mpsqr(GEN x) { return (typ(x)==t_INT) ? sqri(x) : sqrr(x); }
1125 : INLINE GEN
1126 663939 : mpdiv(GEN x, GEN y)
1127 : {
1128 663939 : if (typ(x)==t_INT)
1129 260438 : return (typ(y)==t_INT) ? divii(x,y) : divir(x,y);
1130 403501 : return (typ(y)==t_INT) ? divri(x,y) : divrr(x,y);
1131 : }
1132 :
1133 : /*******************************************************************/
1134 : /* */
1135 : /* Z/nZ, n ULONG */
1136 : /* */
1137 : /*******************************************************************/
1138 : INLINE ulong
1139 439338504 : Fl_double(ulong a, ulong p)
1140 : {
1141 439338504 : ulong res = a << 1;
1142 439338504 : return (res >= p || res < a) ? res - p : res;
1143 : }
1144 : INLINE ulong
1145 88557707 : Fl_triple(ulong a, ulong p)
1146 : {
1147 88557707 : ulong res = a << 1;
1148 88557707 : if (res >= p || res < a) res -= p;
1149 88557707 : res += a;
1150 88557707 : return (res >= p || res < a)? res - p: res;
1151 : }
1152 : INLINE ulong
1153 16894940 : Fl_halve(ulong a, ulong p)
1154 : {
1155 : ulong ap, ap2;
1156 16894940 : if ((a&1UL)==0) return a>>1;
1157 8501052 : ap = a + p; ap2 = ap>>1;
1158 8501052 : return ap>=a ? ap2: (ap2|HIGHBIT);
1159 : }
1160 :
1161 : INLINE ulong
1162 4046460238 : Fl_add(ulong a, ulong b, ulong p)
1163 : {
1164 4046460238 : ulong res = a + b;
1165 4046460238 : return (res >= p || res < a) ? res - p : res;
1166 : }
1167 : INLINE ulong
1168 695374637 : Fl_neg(ulong x, ulong p) { return x ? p - x: 0; }
1169 :
1170 : INLINE ulong
1171 6682540676 : Fl_sub(ulong a, ulong b, ulong p)
1172 : {
1173 6682540676 : ulong res = a - b;
1174 6682540676 : return (res > a) ? res + p: res;
1175 : }
1176 :
1177 : /* centerlift(u mod p) */
1178 : INLINE long
1179 3900221 : Fl_center(ulong u, ulong p, ulong ps2) { return (long) (u > ps2)? u - p: u; }
1180 :
1181 : INLINE ulong
1182 2240664208 : Fl_mul(ulong a, ulong b, ulong p)
1183 : {
1184 : ulong x;
1185 : LOCAL_HIREMAINDER;
1186 2240664208 : x = mulll(a,b);
1187 2240664208 : if (!hiremainder) return x % p;
1188 349512300 : (void)divll(x,p); return hiremainder;
1189 : }
1190 : INLINE ulong
1191 94223330 : Fl_sqr(ulong a, ulong p)
1192 : {
1193 : ulong x;
1194 : LOCAL_HIREMAINDER;
1195 94223330 : x = mulll(a,a);
1196 94223330 : if (!hiremainder) return x % p;
1197 25790798 : (void)divll(x,p); return hiremainder;
1198 : }
1199 : /* don't assume that p is prime: can't special case a = 0 */
1200 : INLINE ulong
1201 32219864 : Fl_div(ulong a, ulong b, ulong p)
1202 32219864 : { return Fl_mul(a, Fl_inv(b, p), p); }
1203 :
1204 : /*******************************************************************/
1205 : /* */
1206 : /* DEFINED FROM EXISTING ONE EXPLOITING COMMUTATIVITY */
1207 : /* */
1208 : /*******************************************************************/
1209 : INLINE GEN
1210 1098741 : addri(GEN x, GEN y) { return addir(y,x); }
1211 : INLINE GEN
1212 145810793 : addis(GEN x, long s) { return addsi(s,x); }
1213 : INLINE GEN
1214 88684645 : addiu(GEN x, ulong s) { return addui(s,x); }
1215 : INLINE GEN
1216 10865024 : addrs(GEN x, long s) { return addsr(s,x); }
1217 :
1218 : INLINE GEN
1219 122351342 : subiu(GEN x, long y) { GEN z = subui(y, x); togglesign(z); return z; }
1220 : INLINE GEN
1221 169856 : subis(GEN x, long y) { return addsi(-y,x); }
1222 : INLINE GEN
1223 13238383 : subrs(GEN x, long y) { return addsr(-y,x); }
1224 :
1225 : INLINE GEN
1226 482365713 : mulis(GEN x, long s) { return mulsi(s,x); }
1227 : INLINE GEN
1228 363797973 : muliu(GEN x, ulong s) { return mului(s,x); }
1229 : INLINE GEN
1230 2680296 : mulru(GEN x, ulong s) { return mulur(s,x); }
1231 : INLINE GEN
1232 34152896 : mulri(GEN x, GEN s) { return mulir(s,x); }
1233 : INLINE GEN
1234 7114837 : mulrs(GEN x, long s) { return mulsr(s,x); }
1235 :
1236 : /*******************************************************************/
1237 : /* */
1238 : /* VALUATION, EXPONENT, SHIFTS */
1239 : /* */
1240 : /*******************************************************************/
1241 : INLINE long
1242 174747242 : vali(GEN x)
1243 : {
1244 : long i;
1245 : GEN xp;
1246 :
1247 174747242 : if (!signe(x)) return -1;
1248 174670449 : xp=int_LSW(x);
1249 179811576 : for (i=0; !*xp; i++) xp=int_nextW(xp);
1250 174670449 : return vals(*xp) + i * BITS_IN_LONG;
1251 : }
1252 :
1253 : /* assume x > 0 */
1254 : INLINE long
1255 678813830 : expu(ulong x) { return (BITS_IN_LONG-1) - (long)bfffo(x); }
1256 :
1257 : INLINE long
1258 1909763929 : expi(GEN x)
1259 : {
1260 1909763929 : const long lx=lgefint(x);
1261 1909763929 : return lx==2? -(long)HIGHEXPOBIT: bit_accuracy(lx)-(long)bfffo(*int_MSW(x))-1;
1262 : }
1263 :
1264 : INLINE GEN
1265 145720723 : shiftr(GEN x, long n)
1266 : {
1267 145720723 : const long e = evalexpo(expo(x)+n);
1268 145716946 : const GEN y = rcopy(x);
1269 :
1270 145705239 : if (e & ~EXPOBITS) pari_err_OVERFLOW("expo()");
1271 145706977 : y[1] = (y[1]&~EXPOBITS) | e; return y;
1272 : }
1273 : INLINE GEN
1274 105739551 : mpshift(GEN x,long s) { return (typ(x)==t_INT)?shifti(x,s):shiftr(x,s); }
1275 :
1276 : /* FIXME: adapt/use mpn_[lr]shift instead */
1277 : /* z2[imin..imax] := z1[imin..imax].f shifted left sh bits
1278 : * (feeding f from the right). Assume sh > 0 */
1279 : INLINE void
1280 6376548391 : shift_left(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh)
1281 : {
1282 6376548391 : GEN sb = z1 + imin, se = z1 + imax, te = z2 + imax;
1283 6376548391 : ulong l, m = BITS_IN_LONG - sh, k = f >> m;
1284 40677767605 : while (se > sb) {
1285 34301219214 : l = *se--;
1286 34301219214 : *te-- = (l << sh) | k;
1287 34301219214 : k = l >> m;
1288 : }
1289 6376548391 : *te = (((ulong)*se) << sh) | k;
1290 6376548391 : }
1291 : /* z2[imin..imax] := f.z1[imin..imax-1] shifted right sh bits
1292 : * (feeding f from the left). Assume sh > 0 */
1293 : INLINE void
1294 4926002503 : shift_right(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh)
1295 : {
1296 4926002503 : GEN sb = z1 + imin, se = z1 + imax, tb = z2 + imin;
1297 4926002503 : ulong k, l = *sb++, m = BITS_IN_LONG - sh;
1298 4926002503 : *tb++ = (l >> sh) | (f << m);
1299 23733765396 : while (sb < se) {
1300 18807762893 : k = l << m;
1301 18807762893 : l = *sb++;
1302 18807762893 : *tb++ = (l >> sh) | k;
1303 : }
1304 4926002503 : }
1305 :
1306 : /* Backward compatibility. Inefficient && unused */
1307 : extern ulong hiremainder;
1308 : INLINE ulong
1309 0 : shiftl(ulong x, ulong y)
1310 0 : { hiremainder = x>>(BITS_IN_LONG-y); return (x<<y); }
1311 :
1312 : INLINE ulong
1313 0 : shiftlr(ulong x, ulong y)
1314 0 : { hiremainder = x<<(BITS_IN_LONG-y); return (x>>y); }
1315 :
1316 : INLINE void
1317 297287726 : shiftr_inplace(GEN z, long d)
1318 : {
1319 297287726 : setexpo(z, expo(z)+d);
1320 297283822 : }
1321 :
1322 : /*******************************************************************/
1323 : /* */
1324 : /* ASSIGNMENT */
1325 : /* */
1326 : /*******************************************************************/
1327 : INLINE void
1328 570664299 : affii(GEN x, GEN y)
1329 : {
1330 570664299 : long lx = lgefint(x);
1331 570664299 : if (lg(y)<lx) pari_err_OVERFLOW("t_INT-->t_INT assignment");
1332 36521879816 : while (--lx) y[lx] = x[lx];
1333 570669037 : }
1334 : INLINE void
1335 4194890 : affsi(long s, GEN x)
1336 : {
1337 4194890 : if (!s) x[1] = evalsigne(0) | evallgefint(2);
1338 : else
1339 : {
1340 3965626 : if (s > 0) { x[1] = evalsigne( 1) | evallgefint(3); x[2] = s; }
1341 1126603 : else { x[1] = evalsigne(-1) | evallgefint(3); x[2] = -s; }
1342 : }
1343 4194890 : }
1344 : INLINE void
1345 44603025 : affui(ulong u, GEN x)
1346 : {
1347 44603025 : if (!u) x[1] = evalsigne(0) | evallgefint(2);
1348 44563839 : else { x[1] = evalsigne(1) | evallgefint(3); x[2] = u; }
1349 44603025 : }
1350 :
1351 : INLINE void
1352 413725414 : affsr(long x, GEN y)
1353 : {
1354 413725414 : long sh, i, ly = lg(y);
1355 :
1356 413725414 : if (!x)
1357 : {
1358 70896 : y[1] = evalexpo(-bit_accuracy(ly));
1359 70896 : return;
1360 : }
1361 413654518 : if (x < 0) {
1362 7973 : x = -x; sh = bfffo(x);
1363 7973 : y[1] = evalsigne(-1) | _evalexpo((BITS_IN_LONG-1)-sh);
1364 : }
1365 : else
1366 : {
1367 413646545 : sh = bfffo(x);
1368 413646545 : y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
1369 : }
1370 4191771210 : y[2] = ((ulong)x)<<sh; for (i=3; i<ly; i++) y[i]=0;
1371 : }
1372 :
1373 : INLINE void
1374 12215314 : affur(ulong x, GEN y)
1375 : {
1376 12215314 : long sh, i, ly = lg(y);
1377 :
1378 12215314 : if (!x)
1379 : {
1380 1299414 : y[1] = evalexpo(-bit_accuracy(ly));
1381 1299414 : return;
1382 : }
1383 10915900 : sh = bfffo(x);
1384 10915900 : y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
1385 43290293 : y[2] = x<<sh; for (i=3; i<ly; i++) y[i] = 0;
1386 : }
1387 :
1388 : INLINE void
1389 424072 : affiz(GEN x, GEN y) { if (typ(y)==t_INT) affii(x,y); else affir(x,y); }
1390 : INLINE void
1391 0 : affsz(long x, GEN y) { if (typ(y)==t_INT) affsi(x,y); else affsr(x,y); }
1392 : INLINE void
1393 658136 : mpaff(GEN x, GEN y) { if (typ(x)==t_INT) affiz(x, y); else affrr(x,y); }
1394 :
1395 : /*******************************************************************/
1396 : /* */
1397 : /* OPERATION + ASSIGNMENT */
1398 : /* */
1399 : /*******************************************************************/
1400 :
1401 0 : INLINE void addiiz(GEN x, GEN y, GEN z)
1402 0 : { pari_sp av = avma; affii(addii(x,y),z); set_avma(av); }
1403 0 : INLINE void addirz(GEN x, GEN y, GEN z)
1404 0 : { pari_sp av = avma; affrr(addir(x,y),z); set_avma(av); }
1405 0 : INLINE void addriz(GEN x, GEN y, GEN z)
1406 0 : { pari_sp av = avma; affrr(addri(x,y),z); set_avma(av); }
1407 1307078 : INLINE void addrrz(GEN x, GEN y, GEN z)
1408 1307078 : { pari_sp av = avma; affrr(addrr(x,y),z); set_avma(av); }
1409 0 : INLINE void addsiz(long s, GEN y, GEN z)
1410 0 : { pari_sp av = avma; affii(addsi(s,y),z); set_avma(av); }
1411 0 : INLINE void addsrz(long s, GEN y, GEN z)
1412 0 : { pari_sp av = avma; affrr(addsr(s,y),z); set_avma(av); }
1413 0 : INLINE void addssz(long s, long y, GEN z)
1414 0 : { pari_sp av = avma; affii(addss(s,y),z); set_avma(av); }
1415 :
1416 0 : INLINE void diviiz(GEN x, GEN y, GEN z)
1417 0 : { pari_sp av = avma; affii(divii(x,y),z); set_avma(av); }
1418 0 : INLINE void divirz(GEN x, GEN y, GEN z)
1419 0 : { pari_sp av = avma; mpaff(divir(x,y),z); set_avma(av); }
1420 0 : INLINE void divisz(GEN x, long y, GEN z)
1421 0 : { pari_sp av = avma; affii(divis(x,y),z); set_avma(av); }
1422 0 : INLINE void divriz(GEN x, GEN y, GEN z)
1423 0 : { pari_sp av = avma; affrr(divri(x,y),z); set_avma(av); }
1424 511 : INLINE void divrrz(GEN x, GEN y, GEN z)
1425 511 : { pari_sp av = avma; affrr(divrr(x,y),z); set_avma(av); }
1426 0 : INLINE void divrsz(GEN y, long s, GEN z)
1427 0 : { pari_sp av = avma; affrr(divrs(y,s),z); set_avma(av); }
1428 0 : INLINE void divsiz(long x, GEN y, GEN z)
1429 0 : { long junk; affsi(sdivsi_rem(x,y,&junk), z); }
1430 0 : INLINE void divsrz(long s, GEN y, GEN z)
1431 0 : { pari_sp av = avma; mpaff(divsr(s,y),z); set_avma(av); }
1432 0 : INLINE void divssz(long x, long y, GEN z)
1433 0 : { affsi(x/y, z); }
1434 :
1435 0 : INLINE void modisz(GEN y, long s, GEN z)
1436 0 : { affsi(smodis(y,s),z); }
1437 0 : INLINE void modsiz(long s, GEN y, GEN z)
1438 0 : { pari_sp av = avma; affii(modsi(s,y),z); set_avma(av); }
1439 0 : INLINE void modssz(long s, long y, GEN z)
1440 0 : { affsi(smodss(s,y),z); }
1441 :
1442 0 : INLINE void mpaddz(GEN x, GEN y, GEN z)
1443 0 : { pari_sp av = avma; mpaff(mpadd(x,y),z); set_avma(av); }
1444 0 : INLINE void mpsubz(GEN x, GEN y, GEN z)
1445 0 : { pari_sp av = avma; mpaff(mpsub(x,y),z); set_avma(av); }
1446 0 : INLINE void mpmulz(GEN x, GEN y, GEN z)
1447 0 : { pari_sp av = avma; mpaff(mpmul(x,y),z); set_avma(av); }
1448 :
1449 0 : INLINE void muliiz(GEN x, GEN y, GEN z)
1450 0 : { pari_sp av = avma; affii(mulii(x,y),z); set_avma(av); }
1451 0 : INLINE void mulirz(GEN x, GEN y, GEN z)
1452 0 : { pari_sp av = avma; mpaff(mulir(x,y),z); set_avma(av); }
1453 0 : INLINE void mulriz(GEN x, GEN y, GEN z)
1454 0 : { pari_sp av = avma; mpaff(mulri(x,y),z); set_avma(av); }
1455 192514 : INLINE void mulrrz(GEN x, GEN y, GEN z)
1456 192514 : { pari_sp av = avma; affrr(mulrr(x,y),z); set_avma(av); }
1457 0 : INLINE void mulsiz(long s, GEN y, GEN z)
1458 0 : { pari_sp av = avma; affii(mulsi(s,y),z); set_avma(av); }
1459 0 : INLINE void mulsrz(long s, GEN y, GEN z)
1460 0 : { pari_sp av = avma; mpaff(mulsr(s,y),z); set_avma(av); }
1461 0 : INLINE void mulssz(long s, long y, GEN z)
1462 0 : { pari_sp av = avma; affii(mulss(s,y),z); set_avma(av); }
1463 :
1464 0 : INLINE void remiiz(GEN x, GEN y, GEN z)
1465 0 : { pari_sp av = avma; affii(remii(x,y),z); set_avma(av); }
1466 0 : INLINE void remisz(GEN y, long s, GEN z)
1467 0 : { pari_sp av = avma; affii(remis(y,s),z); set_avma(av); }
1468 0 : INLINE void remsiz(long s, GEN y, GEN z)
1469 0 : { pari_sp av = avma; affii(remsi(s,y),z); set_avma(av); }
1470 0 : INLINE void remssz(long s, long y, GEN z)
1471 0 : { pari_sp av = avma; affii(remss(s,y),z); set_avma(av); }
1472 :
1473 28 : INLINE void subiiz(GEN x, GEN y, GEN z)
1474 28 : { pari_sp av = avma; affii(subii(x,y),z); set_avma(av); }
1475 0 : INLINE void subirz(GEN x, GEN y, GEN z)
1476 0 : { pari_sp av = avma; affrr(subir(x,y),z); set_avma(av); }
1477 0 : INLINE void subisz(GEN y, long s, GEN z)
1478 0 : { pari_sp av = avma; affii(addsi(-s,y),z); set_avma(av); }
1479 0 : INLINE void subriz(GEN x, GEN y, GEN z)
1480 0 : { pari_sp av = avma; affrr(subri(x,y),z); set_avma(av); }
1481 1296706 : INLINE void subrrz(GEN x, GEN y, GEN z)
1482 1296706 : { pari_sp av = avma; affrr(subrr(x,y),z); set_avma(av); }
1483 0 : INLINE void subrsz(GEN y, long s, GEN z)
1484 0 : { pari_sp av = avma; affrr(addsr(-s,y),z); set_avma(av); }
1485 0 : INLINE void subsiz(long s, GEN y, GEN z)
1486 0 : { pari_sp av = avma; affii(subsi(s,y),z); set_avma(av); }
1487 0 : INLINE void subsrz(long s, GEN y, GEN z)
1488 0 : { pari_sp av = avma; affrr(subsr(s,y),z); set_avma(av); }
1489 0 : INLINE void subssz(long x, long y, GEN z) { addssz(x,-y,z); }
1490 :
1491 : INLINE void
1492 0 : dvmdssz(long x, long y, GEN z, GEN t) {
1493 0 : pari_sp av = avma;
1494 : long r;
1495 0 : affii(divss_rem(x,y, &r), z); set_avma(av); affsi(r,t);
1496 0 : }
1497 : INLINE void
1498 0 : dvmdsiz(long x, GEN y, GEN z, GEN t) {
1499 0 : pari_sp av = avma;
1500 : long r;
1501 0 : affii(divsi_rem(x,y, &r), z); set_avma(av); affsi(r,t);
1502 0 : }
1503 : INLINE void
1504 0 : dvmdisz(GEN x, long y, GEN z, GEN t) {
1505 0 : pari_sp av = avma;
1506 : long r;
1507 0 : affii(divis_rem(x,y, &r),z); set_avma(av); affsi(r,t);
1508 0 : }
1509 : INLINE void
1510 0 : dvmdiiz(GEN x, GEN y, GEN z, GEN t) {
1511 0 : pari_sp av = avma;
1512 : GEN r;
1513 0 : affii(dvmdii(x,y,&r),z); affii(r,t); set_avma(av);
1514 0 : }
|
