1 #line 2 "../src/kernel/none/level1.h"
2 /* $Id: level1.h 7892 2006-04-19 16:18:26Z kb $
3 
4 Copyright (C) 2000  The PARI group.
5 
6 This file is part of the PARI/GP package.
7 
8 PARI/GP is free software; you can redistribute it and/or modify it under the
9 terms of the GNU General Public License as published by the Free Software
10 Foundation. It is distributed in the hope that it will be useful, but WITHOUT
11 ANY WARRANTY WHATSOEVER.
12 
13 Check the License for details. You should have received a copy of it, along
14 with the package; see the file 'COPYING'. If not, write to the Free Software
15 Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */
16 
17 /* This file defines "level 1" kernel functions.
18  * These functions can be inline; if not they are defined externally in
19  * level1.c, which includes this file and never needs to be changed
20  * The following lines are necessary for level0.c and level1.c */
21 
22 #if !defined(INLINE)
23 GEN    mkcol(GEN x);
24 GEN    mkcol2(GEN x, GEN y);
25 GEN    mkcolcopy(GEN x);
26 GEN    mkcomplex(GEN x, GEN y);
27 GEN    mkfrac(GEN x, GEN y);
28 GEN    mkintmod(GEN x, GEN y);
29 GEN    mkpolmod(GEN x, GEN y);
30 GEN    mkintmodu(ulong x, ulong y);
31 GEN    mkmat(GEN x);
32 GEN    mkmat2(GEN x, GEN y);
33 GEN    mkmatcopy(GEN x);
34 GEN    mkrfrac(GEN x, GEN y);
35 GEN    mkvec(GEN x);
36 GEN    mkvec2(GEN x, GEN y);
37 GEN    mkvec2s(long x, long y);
38 GEN    mkvec2copy(GEN x, GEN y);
39 GEN    mkvec3(GEN x, GEN y, GEN z);
40 GEN    mkvec3s(long x, long y, long z);
41 GEN    mkvec4(GEN x, GEN y, GEN z, GEN t);
42 GEN    mkveccopy(GEN x);
43 GEN    mkvecs(long x);
44 GEN    mkvecsmall(long x);
45 GEN    mkvecsmall2(long x,long y);
46 GEN    mkvecsmall3(long x,long y, long z);
47 void   affiz(GEN x, GEN y);
48 void   affsz(long x, GEN y);
49 GEN    addii(GEN x, GEN y);
50 GEN    addir(GEN x, GEN y);
51 GEN    addrr(GEN x, GEN y);
52 GEN    addsi(long x, GEN y);
53 void   affii(GEN x, GEN y);
54 void   affsi(long s, GEN x);
55 void   affsr(long s, GEN x);
56 void   affui(long s, GEN x);
57 void   affur(ulong s, GEN x);
58 GEN    cgetc(long x);
59 GEN    cgetg_copy(long lx, GEN x);
60 GEN    cgetg(long x, long y);
61 GEN    cgeti(long x);
62 GEN    cgetr(long x);
63 int    cmpir(GEN x, GEN y);
64 int    cmpsr(long x, GEN y);
65 GEN    constant_term(GEN x);
66 GEN    ctofp(GEN x, long prec);
67 void   divrrz(GEN x, GEN y, GEN z);
68 GEN    divsi_rem(long x, GEN y, long *rem);
69 void   divsiz(long x, GEN y, GEN z);
70 GEN    divss(long x, long y);
71 GEN    divss_rem(long x, long y, long *rem);
72 void   divssz(long x, long y, GEN z);
73 ulong  Fl_add(ulong a, ulong b, ulong p);
74 long   Fl_center(ulong u, ulong p, ulong ps2);
75 ulong  Fl_div(ulong a, ulong b, ulong p);
76 ulong  Fl_mul(ulong a, ulong b, ulong p);
77 ulong  Fl_neg(ulong x, ulong p);
78 ulong  Fl_sqr(ulong a, ulong p);
79 ulong  Fl_sub(ulong a, ulong b, ulong p);
80 int    dvdii(GEN x, GEN y);
81 int    dvdiiz(GEN x, GEN y, GEN z);
82 int    dvdisz(GEN x, long y, GEN z);
83 int    dvdiuz(GEN x, ulong y, GEN z);
84 void   dvmdiiz(GEN x, GEN y, GEN z, GEN t);
85 GEN    dvmdis(GEN x, long y, GEN *z);
86 void   dvmdisz(GEN x, long y, GEN z, GEN t);
87 GEN    dvmdsi(long x, GEN y, GEN *z);
88 void   dvmdsiz(long x, GEN y, GEN z, GEN t);
89 GEN    dvmdss(long x, long y, GEN *z);
90 void   dvmdssz(long x, long y, GEN z, GEN t);
91 long   evalexpo(long x);
92 long   evallg(long x);
93 long   evalvalp(long x);
94 long   expi(GEN x);
95 GEN    fractor(GEN x, long prec);
96 double gtodouble(GEN x);
97 GEN    icopy_av(GEN x, GEN y);
98 GEN    icopy(GEN x);
99 GEN    init_gen_op(GEN x, long tx, long *lx, long *i);
100 GEN    itor(GEN x, long prec);
101 long   itos(GEN x);
102 long   itos_or_0(GEN x);
103 ulong  itou(GEN x);
104 ulong  itou_or_0(GEN x);
105 GEN    leading_term(GEN x);
106 long   maxss(long x, long y);
107 long   minss(long x, long y);
108 GEN    modis(GEN x, long y);
109 GEN    modsi(long x, GEN y);
110 GEN    modss(long x, long y);
111 GEN    mpabs(GEN x);
112 GEN    mpadd(GEN x, GEN y);
113 void   mpaff(GEN x, GEN y);
114 GEN    mpceil(GEN x);
115 int    mpcmp(GEN x, GEN y);
116 GEN    mpcopy(GEN x);
117 GEN    mpdiv(GEN x, GEN y);
118 GEN    mpfloor(GEN x);
119 GEN    mpmul(GEN x, GEN y);
120 GEN    mpneg(GEN x);
121 GEN    mpround(GEN x);
122 GEN    mpsub(GEN x, GEN y);
123 GEN    mptrunc(GEN x);
124 GEN    new_chunk(size_t x);
125 long   random_bits(long k);
126 GEN    rdivii(GEN x, GEN y, long prec);
127 GEN    rdiviiz(GEN x, GEN y, GEN z);
128 GEN    rdivis(GEN x, long y, long prec);
129 GEN    rdivsi(long x, GEN y, long prec);
130 GEN    rdivss(long x, long y, long prec);
131 GEN    real2n(long n, long prec);
132 GEN    real_1(long prec);
133 GEN    real_m1(long prec);
134 GEN    real_0_bit(long bitprec);
135 GEN    real_0(long prec);
136 void   remiiz(GEN x, GEN y, GEN z);
137 GEN    remis(GEN x, long y);
138 GEN    remsi(long x, GEN y);
139 GEN    remss(long x, long y);
140 GEN    rtor(GEN x, long prec);
141 long   sdivsi(long x, GEN y)
142 long   sdivsi_rem(long x, GEN y, long *rem);
143 long   sdivss_rem(long x, long y, long *rem);
144 void   shift_left2(GEN z2, GEN z1, long min, long M, ulong f, ulong sh, ulong m);
145 void   shift_right2(GEN z2, GEN z1, long min, long M, ulong f, ulong sh, ulong m);
146 ulong  shiftl(ulong x, ulong y);
147 ulong  shiftlr(ulong x, ulong y);
148 GEN    shiftr(GEN x, long n);
149 long   smodis(GEN x, long y);
150 long   smodss(long x, long y);
151 void   stackdummy(pari_sp av, pari_sp ltop);
152 GEN    stoi(long x);
153 GEN    stor(long x, long prec);
154 GEN    subii(GEN x, GEN y);
155 GEN    subir(GEN x, GEN y);
156 GEN    subri(GEN x, GEN y);
157 GEN    subrr(GEN x, GEN y);
158 GEN    subsi(long x, GEN y);
159 ulong  umodui(ulong x, GEN y);
160 GEN    utoi(ulong x);
161 GEN    utoineg(ulong x);
162 GEN    utoipos(ulong x);
163 GEN    utor(ulong s, long prec);
164 long   vali(GEN x);
165 GEN    zerocol(long n);
166 GEN    zeromat(long m, long n);
167 GEN    zeromatcopy(long m, long n);
168 GEN    zeropadic(GEN p, long e);
169 GEN    zeropol(long v);
170 GEN    zeroser(long v, long e);
171 GEN    zerovec(long n);
172 GEN    col_ei(long n, long i);
173 GEN    vec_ei(long n, long i);
174 
175 #else
176 
177 INLINE long
evallg(long x)178 evallg(long x)
179 {
180   if (x & ~LGBITS) pari_err(errlg);
181   return _evallg(x);
182 }
183 
184 INLINE long
evalvalp(long x)185 evalvalp(long x)
186 {
187   const long v = _evalvalp(x);
188   if (v & ~VALPBITS) pari_err(errvalp);
189   return v;
190 }
191 
192 INLINE long
evalexpo(long x)193 evalexpo(long x)
194 {
195   const long v = _evalexpo(x);
196   if (v & ~EXPOBITS) pari_err(errexpo);
197   return v;
198 }
199 
200 INLINE GEN
constant_term(GEN x)201 constant_term(GEN x) { return signe(x)? gel(x,2): gen_0; }
202 INLINE GEN
leading_term(GEN x)203 leading_term(GEN x) { return lg(x) == 2? gen_0: gel(x,lg(x)-1); }
204 
205 /* Inhibit some area gerepile-wise: declare it to be a non recursive
206  * type, of length l. Thus gerepile won't inspect the zone, just copy it.
207  * For the following situation:
208  *   z = cgetg(t,a); av = avma; garbage(); ltop = avma;
209  *   for (i=1; i<HUGE; i++) gel(z,i) = blah();
210  *   stackdummy(av,ltop);
211  * loses (av-ltop) words but save a costly gerepile.
212  */
213 INLINE void
stackdummy(pari_sp av,pari_sp ltop)214 stackdummy(pari_sp av, pari_sp ltop) {
215   long l = ((GEN)av) - ((GEN)ltop);
216   if (l > 0) {
217     GEN z = (GEN)ltop;
218     z[0] = evaltyp(t_VECSMALL) | evallg(l);
219 #ifdef DEBUG
220     { long i; for (i = 1; i < l; i++) z[i] = 0; }
221 #endif
222   }
223 }
224 
225 INLINE GEN
new_chunk(size_t x)226 new_chunk(size_t x) /* x is a number of bytes */
227 {
228   const GEN z = ((GEN) avma) - x;
229   if (x > ((avma-bot)>>TWOPOTBYTES_IN_LONG)) pari_err(errpile);
230 #if defined(_WIN32) || defined(__CYGWIN32__)
231   if (win32ctrlc) dowin32ctrlc();
232 #endif
233   avma = (pari_sp)z;
234 
235 #ifdef MEMSTEP
236   if (DEBUGMEM && memused != DISABLE_MEMUSED) {
237     long d = (long)memused - (long)z;
238     if (d > 4*MEMSTEP || d < -4*MEMSTEP)
239     {
240       memused = (pari_sp)z;
241       fprintferr("...%4.0lf Mbytes used\n",(top-memused)/1048576.);
242     }
243   }
244 #endif
245   return z;
246 }
247 
248 /* cgetg(lg(x), typ(x)), assuming lx = lg(x). Implicit unsetisclone() */
249 INLINE GEN
cgetg_copy(long lx,GEN x)250 cgetg_copy(long lx, GEN x) {
251   GEN y = new_chunk((size_t)lx);
252   y[0] = x[0] & (TYPBITS|LGBITS); return y;
253 }
254 INLINE GEN
init_gen_op(GEN x,long tx,long * lx,long * i)255 init_gen_op(GEN x, long tx, long *lx, long *i) {
256   GEN y;
257   *lx = lg(x); y = cgetg_copy(*lx, x);
258   if (lontyp[tx] == 1) *i = 1; else { y[1] = x[1]; *i = 2; }
259   return y;
260 }
261 
262 INLINE GEN
cgetg(long x,long y)263 cgetg(long x, long y)
264 {
265   const GEN z = new_chunk((size_t)x);
266   z[0] = evaltyp(y) | evallg(x);
267   return z;
268 }
269 
270 INLINE GEN
cgeti(long x)271 cgeti(long x)
272 {
273   const GEN z = new_chunk((size_t)x);
274   z[0] = evaltyp(t_INT) | evallg(x);
275   return z;
276 }
277 
278 INLINE GEN
cgetr(long x)279 cgetr(long x)
280 {
281   const GEN z = new_chunk((size_t)x);
282   z[0] = evaltyp(t_REAL) | evallg(x);
283   return z;
284 }
285 INLINE GEN
cgetc(long l)286 cgetc(long l)
287 {
288   GEN u = cgetg(3,t_COMPLEX); gel(u,1) = cgetr(l); gel(u,2) = cgetr(l);
289   return u;
290 }
291 INLINE GEN
mkintmod(GEN x,GEN y)292 mkintmod(GEN x, GEN y) { GEN v = cgetg(3, t_INTMOD);
293   gel(v,1) = y; gel(v,2) = x; return v; }
294 INLINE GEN
mkpolmod(GEN x,GEN y)295 mkpolmod(GEN x, GEN y) { GEN v = cgetg(3, t_POLMOD);
296   gel(v,1) = y; gel(v,2) = x; return v; }
297 INLINE GEN
mkfrac(GEN x,GEN y)298 mkfrac(GEN x, GEN y) { GEN v = cgetg(3, t_FRAC);
299   gel(v,1) = x; gel(v,2) = y; return v; }
300 INLINE GEN
mkrfrac(GEN x,GEN y)301 mkrfrac(GEN x, GEN y) { GEN v = cgetg(3, t_RFRAC);
302   gel(v,1) = x; gel(v,2) = y; return v; }
303 INLINE GEN
mkcomplex(GEN x,GEN y)304 mkcomplex(GEN x, GEN y) { GEN v = cgetg(3, t_COMPLEX);
305   gel(v,1) = x; gel(v,2) = y; return v; }
306 INLINE GEN
mkvec(GEN x)307 mkvec(GEN x) { GEN v = cgetg(2, t_VEC); gel(v,1) = x; return v; }
308 INLINE GEN
mkvecsmall(long x)309 mkvecsmall(long x) { GEN v = cgetg(2, t_VECSMALL); v[1] = x; return v; }
310 INLINE GEN
mkvecsmall2(long x,long y)311 mkvecsmall2(long x,long y) { GEN v = cgetg(3, t_VECSMALL);
312   v[1]=x; v[2]=y; return v; }
313 INLINE GEN
mkvecsmall3(long x,long y,long z)314 mkvecsmall3(long x,long y, long z) { GEN v = cgetg(4, t_VECSMALL);
315   v[1]=x; v[2]=y; v[3]=z; return v; }
316 INLINE GEN
mkveccopy(GEN x)317 mkveccopy(GEN x) { GEN v = cgetg(2, t_VEC); gel(v,1) = gcopy(x); return v; }
318 INLINE GEN
mkvec2(GEN x,GEN y)319 mkvec2(GEN x, GEN y) {
320   GEN v = cgetg(3,t_VEC); gel(v,1) = x; gel(v,2) = y; return v; }
321 INLINE GEN
mkvec3(GEN x,GEN y,GEN z)322 mkvec3(GEN x, GEN y, GEN z) {
323   GEN v=cgetg(4,t_VEC); gel(v,1) = x; gel(v,2) = y; gel(v,3) = z; return v; }
324 INLINE GEN
mkvec4(GEN x,GEN y,GEN z,GEN t)325 mkvec4(GEN x, GEN y, GEN z, GEN t) {
326   GEN v=cgetg(5,t_VEC); gel(v,1) = x; gel(v,2) = y; gel(v,3) = z; gel(v,4) = t;
327   return v; }
328 INLINE GEN
mkvec2copy(GEN x,GEN y)329 mkvec2copy(GEN x, GEN y) {
330   GEN v = cgetg(3,t_VEC); gel(v,1) = gcopy(x); gel(v,2) = gcopy(y); return v; }
331 INLINE GEN
mkcol(GEN x)332 mkcol(GEN x) { GEN v = cgetg(2, t_COL); gel(v,1) = x; return v; }
333 INLINE GEN
mkcol2(GEN x,GEN y)334 mkcol2(GEN x, GEN y) {
335   GEN v = cgetg(3,t_COL); gel(v,1) = x; gel(v,2) = y; return v; }
336 INLINE GEN
mkcolcopy(GEN x)337 mkcolcopy(GEN x) { GEN v = cgetg(2, t_COL); gel(v,1) = gcopy(x); return v; }
338 INLINE GEN
mkmat(GEN x)339 mkmat(GEN x) { GEN v = cgetg(2, t_MAT); gel(v,1) = x; return v; }
340 INLINE GEN
mkmat2(GEN x,GEN y)341 mkmat2(GEN x, GEN y) { GEN v=cgetg(3,t_MAT); gel(v,1)=x; gel(v,2)=y; return v; }
342 INLINE GEN
mkmatcopy(GEN x)343 mkmatcopy(GEN x) { GEN v = cgetg(2, t_MAT); gel(v,1) = gcopy(x); return v; }
344 
345 /***   ZERO   ***/
346 /* O(p^e) */
347 INLINE GEN
zeropadic(GEN p,long e)348 zeropadic(GEN p, long e)
349 {
350   GEN y = cgetg(5,t_PADIC);
351   gel(y,4) = gen_0;
352   gel(y,3) = gen_1;
353   copyifstack(p,y[2]);
354   y[1] = evalvalp(e) | evalprecp(0);
355   return y;
356 }
357 /* O(pol_x[v]^e) */
358 INLINE GEN
zeroser(long v,long e)359 zeroser(long v, long e)
360 {
361   GEN x = cgetg(2, t_SER);
362   x[1] = evalvalp(e) | evalvarn(v); return x;
363 }
364 /* 0 * pol_x[v] */
365 INLINE GEN
zeropol(long v)366 zeropol(long v)
367 {
368   GEN x = cgetg(2,t_POL);
369   x[1] = evalvarn(v); return x;
370 }
371 /* vector(n) */
372 INLINE GEN
zerocol(long n)373 zerocol(long n)
374 {
375   GEN y = cgetg(n+1,t_COL);
376   long i; for (i=1; i<=n; i++) gel(y,i) = gen_0;
377   return y;
378 }
379 /* vectorv(n) */
380 INLINE GEN
zerovec(long n)381 zerovec(long n)
382 {
383   GEN y = cgetg(n+1,t_VEC);
384   long i; for (i=1; i<=n; i++) gel(y,i) = gen_0;
385   return y;
386 }
387 /* matrix(m, n) */
388 INLINE GEN
zeromat(long m,long n)389 zeromat(long m, long n)
390 {
391   GEN y = cgetg(n+1,t_MAT);
392   GEN v = zerocol(m);
393   long i; for (i=1; i<=n; i++) gel(y,i) = v;
394   return y;
395 }
396 /* matrix(m, n) */
397 INLINE GEN
zeromatcopy(long m,long n)398 zeromatcopy(long m, long n)
399 {
400   GEN y = cgetg(n+1,t_MAT);
401   long i; for (i=1; i<=n; i++) gel(y,i) = zerocol(m);
402   return y;
403 }
404 /* i-th vector in the standard basis */
405 INLINE GEN
col_ei(long n,long i)406 col_ei(long n, long i) { GEN e = zerocol(n); gel(e,i) = gen_1; return e; }
407 INLINE GEN
vec_ei(long n,long i)408 vec_ei(long n, long i) { GEN e = zerovec(n); gel(e,i) = gen_1; return e; }
409 
410 /* cannot do memcpy because sometimes x and y overlap */
411 INLINE GEN
mpcopy(GEN x)412 mpcopy(GEN x)
413 {
414   register long lx = lg(x);
415   const GEN y = cgetg_copy(lx, x);
416   while (--lx > 0) y[lx] = x[lx];
417   return y;
418 }
419 
420 INLINE GEN
icopy(GEN x)421 icopy(GEN x)
422 {
423   register long lx = lgefint(x);
424   const GEN y = cgeti(lx);
425   while (--lx > 0) y[lx] = x[lx];
426   return y;
427 }
428 
429 /* copy integer x as if we had avma = av */
430 INLINE GEN
icopy_av(GEN x,GEN y)431 icopy_av(GEN x, GEN y)
432 {
433   register long lx = lgefint(x);
434   register long ly = lx;
435   y -= lx;
436   while (--lx > 0) y[lx]=x[lx];
437   y[0] = evaltyp(t_INT)|evallg(ly);
438   return y;
439 }
440 
441 INLINE GEN
mpneg(GEN x)442 mpneg(GEN x)
443 {
444   const GEN y=mpcopy(x);
445   setsigne(y,-signe(x)); return y;
446 }
447 
448 INLINE GEN
mpabs(GEN x)449 mpabs(GEN x)
450 {
451   const GEN y=mpcopy(x);
452   if (signe(x)<0) setsigne(y,1);
453   return y;
454 }
455 
456 INLINE long
smodis(GEN x,long y)457 smodis(GEN x, long y)
458 {
459   long rem;
460   const pari_sp av=avma; (void)divis_rem(x,y, &rem); avma=av;
461   return (rem >= 0) ? rem: labs(y) + rem;
462 }
463 
464 INLINE long
smodss(long x,long y)465 smodss(long x, long y)
466 {
467   long rem = x%y;
468   return (rem >= 0)? rem: labs(y) + rem;
469 }
470 
471 /* assume x != 0, return -x as a t_INT */
472 INLINE GEN
utoineg(ulong x)473 utoineg(ulong x)
474 {
475   GEN y = cgeti(3);
476   y[1] = evalsigne(-1)| evallgefint(3); y[2] = x; return y;
477 }
478 /* assume x != 0, return utoi(x) */
479 INLINE GEN
utoipos(ulong x)480 utoipos(ulong x)
481 {
482   GEN y = cgeti(3);
483   y[1] = evalsigne(1)| evallgefint(3); y[2] = x; return y;
484 }
485 INLINE GEN
utoi(ulong x)486 utoi(ulong x) { return x? utoipos(x): gen_0; }
487 INLINE GEN
stoi(long x)488 stoi(long x)
489 {
490   if (!x) return gen_0;
491   return x > 0? utoipos((ulong)x): utoineg((ulong)-x);
492 }
493 INLINE GEN
mkvecs(long x)494 mkvecs(long x) { GEN v = cgetg(2, t_VEC); gel(v,1) = stoi(x); return v; }
495 INLINE GEN
mkvec2s(long x,long y)496 mkvec2s(long x, long y) {
497   GEN v = cgetg(3,t_VEC); gel(v,1) = stoi(x); gel(v,2) = stoi(y); return v; }
498 INLINE GEN
mkvec3s(long x,long y,long z)499 mkvec3s(long x, long y, long z) {
500   GEN v=cgetg(4,t_VEC);
501   gel(v,1)=stoi(x); gel(v,2)=stoi(y); gel(v,3)=stoi(z); return v;
502 }
503 INLINE GEN
mkintmodu(ulong x,ulong y)504 mkintmodu(ulong x, ulong y) {
505   GEN v = cgetg(3,t_INTMOD);
506   gel(v,1) = utoipos(y);
507   gel(v,2) = utoi(x); return v;
508 }
509 
510 INLINE GEN
stosmall(long x)511 stosmall(long x)
512 {
513   if (labs(x) & SMALL_MASK) return stoi(x);
514   return (GEN) (1 | (x<<1));
515 }
516 
517 INLINE long
itos(GEN x)518 itos(GEN x)
519 {
520   const long s = signe(x);
521   long u;
522 
523   if (!s) return 0;
524   u = x[2]; if (lgefint(x) > 3 || u < 0) pari_err(affer2);
525   return (s>0) ? u : -u;
526 }
527 
528 /* as itos, but return 0 if too large. Cf is_bigint */
529 INLINE long
itos_or_0(GEN x)530 itos_or_0(GEN x) {
531   long n;
532   if (lgefint(x) != 3 || (n = x[2]) & HIGHBIT) return 0;
533   return signe(x) > 0? n: -n;
534 }
535 /* as itou, but return 0 if too large. Cf is_bigint */
536 INLINE ulong
itou_or_0(GEN x)537 itou_or_0(GEN x) {
538   if (lgefint(x) != 3) return 0;
539   return (ulong)x[2];
540 }
541 
542 INLINE GEN
modss(long x,long y)543 modss(long x, long y) { return stoi(smodss(x, y)); }
544 
545 INLINE GEN
remss(long x,long y)546 remss(long x, long y) { return stoi(x % y); }
547 
548 INLINE void
affii(GEN x,GEN y)549 affii(GEN x, GEN y)
550 {
551   long lx;
552 
553   if (x==y) return;
554   lx=lgefint(x); if (lg(y)<lx) pari_err(affer3);
555   while (--lx) y[lx]=x[lx];
556 }
557 
558 INLINE void
affsi(long s,GEN x)559 affsi(long s, GEN x)
560 {
561   if (!s) x[1] = evalsigne(0) | evallgefint(2);
562   else
563   {
564     if (s > 0) { x[1] = evalsigne( 1) | evallgefint(3); x[2] =  s; }
565     else       { x[1] = evalsigne(-1) | evallgefint(3); x[2] = -s; }
566   }
567 }
568 
569 INLINE void
affsr(long x,GEN y)570 affsr(long x, GEN y)
571 {
572   long sh, i, ly = lg(y);
573 
574   if (!x)
575   {
576     y[1] = evalexpo(-bit_accuracy(ly));
577     return;
578   }
579   if (x < 0) {
580     x = -x; sh = bfffo(x);
581     y[1] = evalsigne(-1) | _evalexpo((BITS_IN_LONG-1)-sh);
582   }
583   else
584   {
585     sh = bfffo(x);
586     y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
587   }
588   y[2] = x<<sh; for (i=3; i<ly; i++) y[i]=0;
589 }
590 
591 INLINE void
affur(ulong x,GEN y)592 affur(ulong x, GEN y)
593 {
594   long sh, i, ly = lg(y);
595 
596   if (!x)
597   {
598     y[1] = evalexpo(-bit_accuracy(ly));
599     return;
600   }
601   sh = bfffo(x);
602   y[1] = evalsigne(1) | _evalexpo((BITS_IN_LONG-1)-sh);
603   y[2] = x<<sh; for (i=3; i<ly; i++) y[i] = 0;
604 }
605 
606 INLINE void
affiz(GEN x,GEN y)607 affiz(GEN x, GEN y) { if (typ(y)==t_INT) affii(x,y); else affir(x,y); }
608 INLINE void
affsz(long x,GEN y)609 affsz(long x, GEN y) { if (typ(y)==t_INT) affsi(x,y); else affsr(x,y); }
610 INLINE void
mpaff(GEN x,GEN y)611 mpaff(GEN x, GEN y) { if (typ(x)==t_INT) affiz(x, y); else affrr(x,y); }
612 
613 INLINE GEN
real_0_bit(long bitprec)614 real_0_bit(long bitprec) { GEN x=cgetr(2); x[1]=evalexpo(bitprec); return x; }
615 INLINE GEN
real_0(long prec)616 real_0(long prec) { return real_0_bit(-bit_accuracy(prec)); }
617 INLINE GEN
real_1(long prec)618 real_1(long prec) {
619   GEN x = cgetr(prec);
620   long i;
621   x[1] = evalsigne(1) | _evalexpo(0);
622   x[2] = (long)HIGHBIT; for (i=3; i<prec; i++) x[i] = 0;
623   return x;
624 }
625 INLINE GEN
real_m1(long prec)626 real_m1(long prec) {
627   GEN x = cgetr(prec);
628   long i;
629   x[1] = evalsigne(-1) | _evalexpo(0);
630   x[2] = (long)HIGHBIT; for (i=3; i<prec; i++) x[i] = 0;
631   return x;
632 }
633 /* 2.^n */
634 INLINE GEN
real2n(long n,long prec)635 real2n(long n, long prec) { GEN z = real_1(prec); setexpo(z, n); return z; }
636 INLINE GEN
stor(long s,long prec)637 stor(long s, long prec) { GEN z = cgetr(prec); affsr(s,z); return z; }
638 INLINE GEN
utor(ulong s,long prec)639 utor(ulong s, long prec){ GEN z = cgetr(prec); affur(s,z); return z; }
640 INLINE GEN
itor(GEN x,long prec)641 itor(GEN x, long prec) { GEN z = cgetr(prec); affir(x,z); return z; }
642 INLINE GEN
rtor(GEN x,long prec)643 rtor(GEN x, long prec) { GEN z = cgetr(prec); affrr(x,z); return z; }
644 INLINE GEN
ctofp(GEN x,long prec)645 ctofp(GEN x, long prec) { GEN z = cgetg(3,t_COMPLEX);
646   gel(z,1) = gtofp(gel(x,1),prec);
647   gel(z,2) = gtofp(gel(x,2),prec); return z;
648 }
649 
650 INLINE GEN
shiftr(GEN x,long n)651 shiftr(GEN x, long n)
652 {
653   const long e = evalexpo(expo(x)+n);
654   const GEN y = rcopy(x);
655 
656   if (e & ~EXPOBITS) pari_err(talker,"overflow in real shift");
657   y[1] = (y[1]&~EXPOBITS) | e; return y;
658 }
659 
660 INLINE int
cmpir(GEN x,GEN y)661 cmpir(GEN x, GEN y)
662 {
663   pari_sp av;
664   GEN z;
665 
666   if (!signe(x)) return -signe(y);
667   if (!signe(y)) return  signe(x);
668   av=avma; z = itor(x, lg(y)); avma=av;
669   return cmprr(z,y); /* cmprr does no memory adjustment */
670 }
671 
672 INLINE int
cmpsr(long x,GEN y)673 cmpsr(long x, GEN y)
674 {
675   pari_sp av;
676   GEN z;
677 
678   if (!x) return -signe(y);
679   av=avma; z = stor(x, 3); avma=av;
680   return cmprr(z,y);
681 }
682 
683 INLINE long
maxss(long x,long y)684 maxss(long x, long y) { return x>y?x:y; }
685 INLINE long
minss(long x,long y)686 minss(long x, long y) { return x<y?x:y; }
687 
688 INLINE GEN
subii(GEN x,GEN y)689 subii(GEN x, GEN y)
690 {
691   if (x==y) return gen_0; /* frequent with x = y = gen_0 */
692   return addii_sign(x, signe(x), y, -signe(y));
693 }
694 INLINE GEN
addii(GEN x,GEN y)695 addii(GEN x, GEN y) { return addii_sign(x, signe(x), y, signe(y)); }
696 INLINE GEN
addrr(GEN x,GEN y)697 addrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, signe(y)); }
698 INLINE GEN
subrr(GEN x,GEN y)699 subrr(GEN x, GEN y) { return addrr_sign(x, signe(x), y, -signe(y)); }
700 INLINE GEN
addir(GEN x,GEN y)701 addir(GEN x, GEN y) { return addir_sign(x, signe(x), y, signe(y)); }
702 INLINE GEN
subir(GEN x,GEN y)703 subir(GEN x, GEN y) { return addir_sign(x, signe(x), y, -signe(y)); }
704 INLINE GEN
subri(GEN x,GEN y)705 subri(GEN x, GEN y) { return addir_sign(y, -signe(y), x, signe(x)); }
706 INLINE GEN
addsi(long x,GEN y)707 addsi(long x, GEN y) { return addsi_sign(x, y, signe(y)); }
708 INLINE GEN
subsi(long x,GEN y)709 subsi(long x, GEN y) { return addsi_sign(x, y, -signe(y)); }
710 
711 INLINE long
vali(GEN x)712 vali(GEN x)
713 {
714   long i;
715   GEN xp;
716 
717   if (!signe(x)) return -1;
718   xp=int_LSW(x);
719   for (i=0; !*xp; i++) xp=int_nextW(xp);
720   return (i<<TWOPOTBITS_IN_LONG) + vals(*xp);
721 }
722 
723 INLINE GEN
divss(long x,long y)724 divss(long x, long y) { return stoi(x / y); }
725 
726 INLINE long
sdivss_rem(long x,long y,long * rem)727 sdivss_rem(long x, long y, long *rem)
728 {
729   long q;
730   LOCAL_HIREMAINDER;
731   if (!y) pari_err(gdiver);
732   hiremainder = 0; q = divll((ulong)labs(x),(ulong)labs(y));
733   if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
734   if (y < 0) q = -q;
735   *rem = hiremainder; return q;
736 }
737 
738 INLINE GEN
divss_rem(long x,long y,long * rem)739 divss_rem(long x, long y, long *rem) { return stoi(sdivss_rem(x,y,rem)); }
740 
741 INLINE GEN
dvmdss(long x,long y,GEN * z)742 dvmdss(long x, long y, GEN *z)
743 {
744   long rem;
745   const GEN q = divss_rem(x,y, &rem);
746   *z = stoi(rem); return q;
747 }
748 
749 INLINE long
sdivsi_rem(long x,GEN y,long * rem)750 sdivsi_rem(long x, GEN y, long *rem)
751 {
752   long q, s = signe(y);
753   LOCAL_HIREMAINDER;
754 
755   if (!s) pari_err(gdiver);
756   if (!x || lgefint(y)>3 || ((long)y[2]) < 0) { *rem = x; return 0; }
757   hiremainder=0; q = (long)divll(labs(x), (ulong)y[2]);
758   if (x < 0) { hiremainder = -((long)hiremainder); q = -q; }
759   if (s < 0) q = -q;
760   *rem = hiremainder; return q;
761 }
762 
763 INLINE long
sdivsi(long x,GEN y)764 sdivsi(long x, GEN y)
765 {
766   long q, s = signe(y);
767 
768   if (!s) pari_err(gdiver);
769   if (!x || lgefint(y)>3 || ((long)y[2]) < 0) return 0;
770   q = labs(x) / y[2];
771   if (x < 0) q = -q;
772   if (s < 0) q = -q;
773   return q;
774 }
775 
776 INLINE GEN
modsi(long x,GEN y)777 modsi(long x, GEN y) {
778   long r;
779   (void)sdivsi_rem(x, y, &r);
780   return (r >= 0)? stoi(r): addsi_sign(r, y, 1);
781 }
782 
783 INLINE GEN
divsi_rem(long s,GEN y,long * rem)784 divsi_rem(long s, GEN y, long *rem) { return stoi(sdivsi_rem(s,y,rem)); }
785 
786 INLINE GEN
dvmdsi(long x,GEN y,GEN * z)787 dvmdsi(long x, GEN y, GEN *z)
788 {
789   long rem;
790   const GEN q = divsi_rem(x,y, &rem);
791   *z = stoi(rem); return q;
792 }
793 
794 INLINE GEN
dvmdis(GEN x,long y,GEN * z)795 dvmdis(GEN x, long y, GEN *z)
796 {
797   long rem;
798   const GEN q = divis_rem(x,y, &rem);
799   *z = stoi(rem); return q;
800 }
801 
802 INLINE void
dvmdssz(long x,long y,GEN z,GEN t)803 dvmdssz(long x, long y, GEN z, GEN t)
804 {
805   long rem;
806   const pari_sp av=avma; affiz(divss_rem(x,y, &rem), z); avma=av;
807   affsi(rem,t);
808 }
809 
810 INLINE void
dvmdsiz(long x,GEN y,GEN z,GEN t)811 dvmdsiz(long x, GEN y, GEN z, GEN t)
812 {
813   long rem;
814   const pari_sp av = avma; affiz(divsi_rem(x,y, &rem), z); avma = av;
815   affsi(rem,t);
816 }
817 
818 INLINE void
dvmdisz(GEN x,long y,GEN z,GEN t)819 dvmdisz(GEN x, long y, GEN z, GEN t)
820 {
821   long rem;
822   const pari_sp av=avma; affiz(divis_rem(x,y, &rem),z); avma=av;
823   affsz(rem,t);
824 }
825 
826 INLINE void
dvmdiiz(GEN x,GEN y,GEN z,GEN t)827 dvmdiiz(GEN x, GEN y, GEN z, GEN t)
828 {
829   const pari_sp av=avma;
830   GEN p;
831 
832   affiz(dvmdii(x,y,&p),z); affiz(p,t); avma=av;
833 }
834 
835 INLINE GEN
modis(GEN x,long y)836 modis(GEN x, long y)
837 {
838   return stoi(smodis(x,y));
839 }
840 
841 INLINE ulong
umodui(ulong x,GEN y)842 umodui(ulong x, GEN y)
843 {
844   LOCAL_HIREMAINDER;
845   if (!signe(y)) pari_err(gdiver);
846   if (!x || lgefint(y) > 3) return x;
847   hiremainder = 0; (void)divll(x, y[2]); return hiremainder;
848 }
849 
850 INLINE GEN
remsi(long x,GEN y)851 remsi(long x, GEN y)
852 {
853   long rem;
854   const pari_sp av=avma; (void)divsi_rem(x,y, &rem); avma=av;
855   return stoi(rem);
856 }
857 
858 INLINE GEN
remis(GEN x,long y)859 remis(GEN x, long y)
860 {
861   long rem;
862   const pari_sp av=avma; (void)divis_rem(x,y, &rem); avma=av;
863   return stoi(rem);
864 }
865 
866 INLINE GEN
rdivis(GEN x,long y,long prec)867 rdivis(GEN x, long y, long prec)
868 {
869   GEN z = cgetr(prec);
870   const pari_sp av = avma;
871   affrr(divrs(itor(x,prec), y),z);
872   avma = av; return z;
873 }
874 
875 INLINE void
divsiz(long x,GEN y,GEN z)876 divsiz(long x, GEN y, GEN z)
877 {
878   long junk;
879   affsi(sdivsi_rem(x,y,&junk), z);
880 }
881 
882 INLINE void
divssz(long x,long y,GEN z)883 divssz(long x, long y, GEN z) { affsi(x/y, z); }
884 
885 INLINE GEN
rdivsi(long x,GEN y,long prec)886 rdivsi(long x, GEN y, long prec)
887 {
888   GEN z = cgetr(prec);
889   const pari_sp av = avma;
890   affrr(divsr(x, itor(y,prec)), z);
891   avma = av; return z;
892 }
893 
894 INLINE GEN
rdivss(long x,long y,long prec)895 rdivss(long x, long y, long prec)
896 {
897   GEN z = cgetr(prec);
898   const pari_sp av = avma;
899   affrr(divrs(stor(x, prec), y), z);
900   avma = av; return z;
901 }
902 
903 INLINE GEN
rdiviiz(GEN x,GEN y,GEN z)904 rdiviiz(GEN x, GEN y, GEN z)
905 {
906   const pari_sp av = avma;
907   long prec = lg(z);
908   affir(x, z);
909   if (!is_bigint(y)) {
910     affrr(divrs(z, y[2]), z);
911     if (signe(y) < 0) setsigne(z, -signe(z));
912   }
913   else
914     affrr(divrr(z, itor(y,prec)), z);
915   avma = av; return z;
916 }
917 
918 INLINE GEN
rdivii(GEN x,GEN y,long prec)919 rdivii(GEN x, GEN y, long prec) { return rdiviiz(x, y, cgetr(prec)); }
920 INLINE GEN
fractor(GEN x,long prec)921 fractor(GEN x, long prec) { return rdivii(gel(x,1), gel(x,2), prec); }
922 
923 INLINE void
divrrz(GEN x,GEN y,GEN z)924 divrrz(GEN x, GEN y, GEN z)
925 {
926   const pari_sp av=avma;
927   affrr(divrr(x,y),z); avma=av;
928 }
929 
930 INLINE void
remiiz(GEN x,GEN y,GEN z)931 remiiz(GEN x, GEN y, GEN z)
932 {
933   const pari_sp av=avma;
934   affii(remii(x,y),z); avma=av;
935 }
936 
937 INLINE int
dvdii(GEN x,GEN y)938 dvdii(GEN x, GEN y)
939 {
940   const pari_sp av=avma;
941   const GEN p1=remii(x,y);
942   avma=av; return p1 == gen_0;
943 }
944 
945 INLINE int
mpcmp(GEN x,GEN y)946 mpcmp(GEN x, GEN y)
947 {
948   if (typ(x)==t_INT)
949     return (typ(y)==t_INT) ? cmpii(x,y) : cmpir(x,y);
950   return (typ(y)==t_INT) ? -cmpir(y,x) : cmprr(x,y);
951 }
952 
953 INLINE GEN
mptrunc(GEN x)954 mptrunc(GEN x) { return typ(x)==t_INT? icopy(x): truncr(x); }
955 INLINE GEN
mpfloor(GEN x)956 mpfloor(GEN x) { return typ(x)==t_INT? icopy(x): floorr(x); }
957 INLINE GEN
mpceil(GEN x)958 mpceil(GEN x) { return typ(x)==t_INT? icopy(x): ceilr(x); }
959 INLINE GEN
mpround(GEN x)960 mpround(GEN x) { return typ(x) == t_INT? icopy(x): roundr(x); }
961 
962 INLINE GEN
mpadd(GEN x,GEN y)963 mpadd(GEN x, GEN y)
964 {
965   if (typ(x)==t_INT)
966     return (typ(y)==t_INT) ? addii(x,y) : addir(x,y);
967   return (typ(y)==t_INT) ? addir(y,x) : addrr(x,y);
968 }
969 
970 INLINE GEN
mpsub(GEN x,GEN y)971 mpsub(GEN x, GEN y)
972 {
973   if (typ(x)==t_INT)
974     return (typ(y)==t_INT) ? subii(x,y) : subir(x,y);
975   return (typ(y)==t_INT) ? subri(x,y) : subrr(x,y);
976 }
977 
978 INLINE GEN
mpmul(GEN x,GEN y)979 mpmul(GEN x, GEN y)
980 {
981   if (typ(x)==t_INT)
982     return (typ(y)==t_INT) ? mulii(x,y) : mulir(x,y);
983   return (typ(y)==t_INT) ? mulir(y,x) : mulrr(x,y);
984 }
985 
986 INLINE GEN
mpdiv(GEN x,GEN y)987 mpdiv(GEN x, GEN y)
988 {
989   if (typ(x)==t_INT)
990     return (typ(y)==t_INT) ? divii(x,y) : divir(x,y);
991   return (typ(y)==t_INT) ? divri(x,y) : divrr(x,y);
992 }
993 
994 INLINE int
dvdiiz(GEN x,GEN y,GEN z)995 dvdiiz(GEN x, GEN y, GEN z)
996 {
997   const pari_sp av=avma;
998   GEN p2;
999   const GEN p1=dvmdii(x,y,&p2);
1000 
1001   if (signe(p2)) { avma=av; return 0; }
1002   affii(p1,z); avma=av; return 1;
1003 }
1004 
1005 /* assume 0 <= k < 32. Return random 0 <= x < (1<<k) */
1006 INLINE long
random_bits(long k)1007 random_bits(long k) { return pari_rand31() >> (31 - k); }
1008 
1009 INLINE ulong
itou(GEN x)1010 itou(GEN x)
1011 {
1012   const long s = signe(x);
1013 
1014   if (!s) return 0;
1015   if (lgefint(x) > 3) pari_err(affer2);
1016   return x[2];
1017 }
1018 
1019 INLINE void
affui(ulong u,GEN x)1020 affui(ulong u, GEN x)
1021 {
1022   if (!u) x[1] = evalsigne(0) | evallgefint(2);
1023   else  { x[1] = evalsigne(1) | evallgefint(3); x[2] = u; }
1024 }
1025 
1026 INLINE int
dvdisz(GEN x,long y,GEN z)1027 dvdisz(GEN x, long y, GEN z)
1028 {
1029   const pari_sp av = avma;
1030   long rem;
1031   GEN p1 = divis_rem(x,y, &rem);
1032   avma = av; if (rem) return 0;
1033   affii(p1,z); return 1;
1034 }
1035 
1036 INLINE int
dvdiuz(GEN x,ulong y,GEN z)1037 dvdiuz(GEN x, ulong y, GEN z)
1038 {
1039   const pari_sp av = avma;
1040   ulong rem;
1041   GEN p1 = diviu_rem(x,y, &rem);
1042   avma = av; if (rem) return 0;
1043   affii(p1,z); return 1;
1044 }
1045 
1046 INLINE double
gtodouble(GEN x)1047 gtodouble(GEN x)
1048 {
1049   static long reel4[4]={ evaltyp(t_REAL) | _evallg(4),0,0,0 };
1050 
1051   if (typ(x)==t_REAL) return rtodbl(x);
1052   gaffect(x,(GEN)reel4); return rtodbl((GEN)reel4);
1053 }
1054 
1055 /* same as Fl_add, assume p <= 2^(BIL - 1), so that overflow can't occur */
1056 INLINE ulong
Fl_add_noofl(ulong a,ulong b,ulong p)1057 Fl_add_noofl(ulong a, ulong b, ulong p)
1058 {
1059   ulong res = a + b;
1060   return (res >= p) ? res - p : res;
1061 }
1062 INLINE ulong
Fl_add(ulong a,ulong b,ulong p)1063 Fl_add(ulong a, ulong b, ulong p)
1064 {
1065   ulong res = a + b;
1066   return (res >= p || res < a) ? res - p : res;
1067 }
1068 INLINE ulong
Fl_neg(ulong x,ulong p)1069 Fl_neg(ulong x, ulong p) { return x ? p - x: 0; }
1070 
1071 INLINE ulong
Fl_sub(ulong a,ulong b,ulong p)1072 Fl_sub(ulong a, ulong b, ulong p)
1073 {
1074   ulong res = a - b;
1075   return (res > a) ? res + p: res;
1076 }
1077 
1078 /* centerlift(u mod p) */
1079 INLINE long
Fl_center(ulong u,ulong p,ulong ps2)1080 Fl_center(ulong u, ulong p, ulong ps2) { return (long) (u > ps2)? u - p: u; }
1081 
1082 INLINE ulong
Fl_mul(ulong a,ulong b,ulong p)1083 Fl_mul(ulong a, ulong b, ulong p)
1084 {
1085   LOCAL_HIREMAINDER;
1086   {
1087     register ulong x = mulll(a,b);
1088     (void)divll(x,p);
1089   }
1090   return hiremainder;
1091 }
1092 INLINE ulong
Fl_sqr(ulong a,ulong p)1093 Fl_sqr(ulong a, ulong p)
1094 {
1095   LOCAL_HIREMAINDER;
1096   {
1097     register ulong x = mulll(a,a);
1098     (void)divll(x,p);
1099   }
1100   return hiremainder;
1101 }
1102 INLINE ulong
Fl_div(ulong a,ulong b,ulong p)1103 Fl_div(ulong a, ulong b, ulong p)
1104 {
1105   return Fl_mul(a, Fl_inv(b, p), p);
1106 }
1107 
1108 INLINE long
expi(GEN x)1109 expi(GEN x)
1110 {
1111   const long lx=lgefint(x);
1112   return lx==2? -(long)HIGHEXPOBIT: bit_accuracy(lx)-(long)bfffo(*int_MSW(x))-1;
1113 }
1114 
1115 /* z2[imin..imax] := z1[imin..imax].f shifted left sh bits
1116  * (feeding f from the right). Assume sh > 0 */
1117 INLINE void
shift_left2(GEN z2,GEN z1,long imin,long imax,ulong f,ulong sh,ulong m)1118 shift_left2(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh, ulong m)
1119 {
1120   GEN sb = z1 + imin, se = z1 + imax, te = z2 + imax;
1121   ulong l, k = f >> m;
1122   while (se > sb) {
1123     l     = *se--;
1124     *te-- = (l << sh) | k;
1125     k     = l >> m;
1126   }
1127   *te = (*se << sh) | k;
1128 }
1129 /* z2[imin..imax] := f.z1[imin..imax-1] shifted right sh bits
1130  * (feeding f from the left). Assume sh > 0 */
1131 INLINE void
shift_right2(GEN z2,GEN z1,long imin,long imax,ulong f,ulong sh,ulong m)1132 shift_right2(GEN z2, GEN z1, long imin, long imax, ulong f, ulong sh, ulong m)
1133 {
1134   GEN sb = z1 + imin, se = z1 + imax, tb = z2 + imin;
1135   ulong k, l = *sb++;
1136   *tb++ = (l >> sh) | (f << m);
1137   while (sb < se) {
1138     k     = l << m;
1139     l     = *sb++;
1140     *tb++ = (l >> sh) | k;
1141   }
1142 }
1143 
1144 /* Backward compatibility. Inefficient && unused */
1145 extern ulong hiremainder;
1146 INLINE ulong
shiftl(ulong x,ulong y)1147 shiftl(ulong x, ulong y)
1148 {
1149   hiremainder=x>>(BITS_IN_LONG-y);
1150   return (x<<y);
1151 }
1152 INLINE ulong
shiftlr(ulong x,ulong y)1153 shiftlr(ulong x, ulong y)
1154 {
1155   hiremainder=x<<(BITS_IN_LONG-y);
1156   return (x>>y);
1157 }
1158 #endif
1159