![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to gen2.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / pari / src / basemath|
Return to gen2.c CVS log | Up to [local] / OpenXM_contrib / pari / src / basemath |
1.1 maekawa 1: /********************************************************************/
2: /** **/
3: /** GENERIC OPERATIONS **/
4: /** (second part) **/
5: /** **/
6: /********************************************************************/
7: /* $Id: gen2.c,v 1.3 1999/09/23 17:50:56 karim Exp $ */
8: #include "pari.h"
9:
10: /*******************************************************************/
11: /* */
12: /* OPERATIONS BY VALUE */
13: /* f is a pointer to the function called. */
14: /* result is gaffect-ed to last parameter */
15: /* */
16: /*******************************************************************/
17:
18: void
19: gop1z(GEN (*f)(GEN), GEN x, GEN y)
20: {
21: long av=avma; gaffect(f(x),y); avma=av;
22: }
23:
24: void
25: gop2z(GEN (*f)(GEN, GEN), GEN x, GEN y, GEN z)
26: {
27: long av=avma; gaffect(f(x,y),z); avma=av;
28: }
29:
30: void
31: gops2gsz(GEN (*f)(GEN, long), GEN x, long s, GEN z)
32: {
33: long av=avma; gaffect(f(x,s),z); avma=av;
34: }
35:
36: void
37: gops2sgz(GEN (*f)(long, GEN), long s, GEN y, GEN z)
38: {
39: long av=avma; gaffect(f(s,y),z); avma=av;
40: }
41:
42: void
43: gops2ssz(GEN (*f)(long, long), long s, long y, GEN z)
44: {
45: long av=avma; gaffect(f(s,y),z); avma=av;
46: }
47:
48: /*******************************************************************/
49: /* */
50: /* OPERATIONS USING SMALL INTEGERS */
51: /* */
52: /*******************************************************************/
53:
54: /* small int prototype */
55: static long court_p[] = { evaltyp(t_INT) | m_evallg(3),0,0 };
56:
57: GEN
58: gopsg2(GEN (*f)(GEN, GEN), long s, GEN y)
59: {
60: affsi(s,court_p); return f(court_p,y);
61: }
62:
63: GEN
64: gopgs2(GEN (*f)(GEN, GEN), GEN y, long s)
65: {
66: affsi(s,court_p); return f(y,court_p);
67: }
68:
69: long
70: opgs2(int (*f)(GEN, GEN), GEN y, long s)
71: {
72: affsi(s,court_p); return f(y,court_p);
73: }
74:
75: void
76: gopsg2z(GEN (*f)(GEN, GEN), long s, GEN y, GEN z)
77: {
78: long av=avma;
79: affsi(s,court_p); gaffect(f(court_p,y),z); avma=av;
80: }
81:
82: void
83: gopgs2z(GEN (*f)(GEN, GEN), GEN y, long s, GEN z)
84: {
85: long av=avma;
86: affsi(s,court_p); gaffect(f(y,court_p),z); avma=av;
87: }
88:
89: /*******************************************************************/
90: /* */
91: /* CREATION OF A P-ADIC GEN */
92: /* */
93: /*******************************************************************/
94:
95: GEN
96: cgetp(GEN x)
97: {
98: GEN y = cgetg(5,t_PADIC);
99: y[1] = evalprecp(precp(x)) | evalvalp(0);
100: icopyifstack(x[2], y[2]);
101: y[3] = licopy((GEN)x[3]);
102: y[4] = lgeti(lgefint(x[3])); return y;
103: }
104:
105: /* y[4] not filled */
106: GEN
107: cgetp2(GEN x, long v)
108: {
109: GEN y = cgetg(5,t_PADIC);
110: y[1] = evalprecp(precp(x)) | evalvalp(v);
111: icopyifstack(x[2], y[2]);
112: y[3] = licopy((GEN)x[3]); return y;
113: }
114:
115: /*******************************************************************/
116: /* */
117: /* CLONING & COPY */
118: /* Replicate an existing GEN */
119: /* */
120: /*******************************************************************/
121: /* lontyp = 0 means non recursive type
122: * otherwise:
123: * lontyp = number of codewords
124: * if not in stack, we don't copy the words in [lontyp,lontyp2[
125: */
126: const long lontyp[] = { 0,0,0,1,1,1,1,2,1,1, 2,2,0,1,1,1,1,1,1,1, 2,0,0 };
127: static long lontyp2[] = { 0,0,0,2,1,1,1,3,2,2, 2,2,0,1,1,1,1,1,1,1, 2,0,0 };
128:
129: /* can't do a memcpy there: avma and x may overlap. memmove is slower */
130: GEN
131: gcopy(GEN x)
132: {
133: long tx=typ(x),lx,i;
134: GEN y;
135:
136: if (tx == t_SMALL) return x;
137: lx = lg(x); y=new_chunk(lx);
138: if (! is_recursive_t(tx))
139: for (i=lx-1; i>=0; i--) y[i]=x[i];
140: else
141: {
142: if (tx==t_POL || tx==t_LIST) lx=lgef(x);
143: for (i=0; i<lontyp[tx]; i++) y[i]=x[i];
144: for ( ; i<lontyp2[tx]; i++) copyifstack(x[i],y[i]);
145: for ( ; i<lx; i++) y[i]=lcopy((GEN)x[i]);
146: }
147: /* unsetisclone(y); useless because gunclone checks isonstack */
148: return y;
149: }
150:
151: GEN
152: gcopy_i(GEN x, long lx)
153: {
154: long tx=typ(x),i;
155: GEN y;
156:
157: if (tx == t_SMALL) return x;
158: y=cgetg(lx,tx);
159: if (! is_recursive_t(tx))
160: for (i=lx-1; i>0; i--) y[i]=x[i];
161: else
162: {
163: for (i=1; i<lontyp[tx]; i++) y[i]=x[i];
164: for ( ; i<lontyp2[tx]; i++) copyifstack(x[i],y[i]);
165: for ( ; i<lx; i++) y[i]=lcopy((GEN)x[i]);
166: }
167: return y;
168: }
169:
170: GEN
171: forcecopy(GEN x)
172: {
173: long tx=typ(x),lx,i;
174: GEN y;
175:
176: if (tx == t_SMALL) return x;
177: lx=lg(x); y=new_chunk(lx);
178: if (! is_recursive_t(tx))
179: for (i=lx-1; i>=0; i--) y[i]=x[i];
180: else
181: {
182: if (tx==t_POL || tx==t_LIST) lx=lgef(x);
183: for (i=0; i<lontyp[tx]; i++) y[i]=x[i];
184: for ( ; i<lx; i++) y[i]=(long)forcecopy((GEN)x[i]);
185: }
186: unsetisclone(y); return y;
187: }
188:
189: GEN
190: dummycopy(GEN x)
191: {
192: long tx=typ(x), lx=lg(x),i;
193: GEN y=new_chunk(lx);
194:
195: switch(tx)
196: {
197: case t_POLMOD:
198: y[1]=x[1]; y[2]=(long) dummycopy((GEN)x[2]);
199: break;
200: case t_MAT:
201: for (i=lx-1;i;i--) y[i]=(long)dummycopy((GEN)x[i]);
202: break;
203: default:
204: for (i=lx-1;i;i--) y[i]=x[i];
205: }
206: y[0]=x[0]; return y;
207: }
208:
209: /* assume x is a t_POL. FOR INTERNAL USE! */
210: GEN
211: dummyclone(GEN x)
212: {
213: long lx = lgef(x);
214: GEN z = (GEN)gpmalloc(lx*sizeof(long));
215: while (--lx >= 0) z[lx] = x[lx];
216: return z;
217: }
218:
219: long
220: taille(GEN x)
221: {
222: long i,n,lx=lg(x),tx=typ(x);
223:
224: n = lx;
225: if (is_recursive_t(tx))
226: {
227: if (tx==t_POL || tx==t_LIST) lx = lgef(x);
228: for (i=lontyp[tx]; i<lx; i++) n += taille((GEN)x[i]);
229: }
230: return n;
231: }
232:
233: long
234: glength(GEN x)
235: {
236: switch(typ(x))
237: {
238: case t_INT: return lgefint(x)-2;
239: case t_POL: case t_LIST: return lgef(x)-2;
240: case t_REAL: return signe(x)? lg(x)-2: 0;
241: case t_STR: return strlen(GSTR(x));
242: }
243: return lg(x)-lontyp[typ(x)];
244: }
245:
246: GEN
247: matsize(GEN x)
248: {
249: GEN y=cgetg(3,t_VEC);
250:
251: switch(typ(x))
252: {
253: case t_VEC:
254: y[1]=un; y[2]=lstoi(lg(x)-1); break;
255: case t_COL:
256: y[1]=lstoi(lg(x)-1); y[2]=un; break;
257: case t_MAT:
258: y[2]=lstoi(lg(x)-1);
259: y[1]=(lg(x)==1)? zero: lstoi(lg(x[1])-1); break;
260: default: err(typeer,"matsize");
261: }
262: return y;
263: }
264:
265: long
266: taille2(GEN x) { return taille(x)<<TWOPOTBYTES_IN_LONG; }
267:
268: GEN
269: brutcopy(GEN x, GEN y)
270: {
271: long i,lx,tx=typ(x);
272: GEN z;
273:
274: if (! is_recursive_t(tx))
275: {
276: lx = (tx==t_INT)? lgefint(x): lg(x);
277: for (i=0; i<lx; i++) y[i] = x[i];
278: }
279: else
280: {
281: lx=lg(x); z = y + lx;
282: if (tx==t_POL || tx==t_LIST) lx = lgef(x);
283: for (i=0; i<lontyp[tx]; i++) y[i] = x[i];
284: for ( ; i<lx; i++)
285: {
286: y[i] = (long)brutcopy((GEN)x[i], z);
287: z += taille((GEN)x[i]);
288: }
289: }
290: unsetisclone(y); return y;
291: }
292:
293: GEN
294: gclone(GEN x)
295: {
296: x = brutcopy(x, newbloc(taille(x)));
297: setisclone(x); return x;
298: }
299:
300: /*******************************************************************/
301: /* */
302: /* GREFFE */
303: /* Greffe d'une serie sur un polynome */
304: /* */
305: /*******************************************************************/
306:
307: GEN
308: greffe(GEN x, long l, long use_stack)
309: {
310: long i,e,k,vx;
311: GEN y;
312:
313: if (typ(x)!=t_POL) err(notpoler,"greffe");
314: if (use_stack) y = cgetg(l,t_SER);
315: else
316: {
317: y = (GEN) gpmalloc(l*sizeof(long));
318: y[0] = evaltyp(t_SER) | evallg(l);
319: }
320: if (gcmp0(x))
321: {
322: y[1]=evalvalp(l-2) | evalvarn(varn(x));
323: i=2; while(i<l) { y[i]=x[2]; i++; }
324: return y;
325: }
326: vx=varn(x); e=gval(x,vx);
327: y[1]=evalsigne(1) | evalvalp(e) | evalvarn(vx);
328: k=lgef(x)-e-1; i=l-1;
329: if (k<l)
330: while (i>k) { y[i]=zero; i--; }
331: while (i>=2) { y[i]=x[i+e]; i--; }
332: return y;
333: }
334:
335: /*******************************************************************/
336: /* */
337: /* CONVERSION GEN --> long */
338: /* */
339: /*******************************************************************/
340:
341: long
342: gtolong(GEN x)
343: {
344: long y,tx=typ(x),av=avma;
345:
346: switch(tx)
347: {
348: case t_INT:
349: return itos(x);
350: case t_REAL: case t_FRAC: case t_FRACN:
351: y=itos(ground(x)); avma=av; return y;
352: case t_COMPLEX:
353: if (gcmp0((GEN)x[2])) return gtolong((GEN)x[1]); break;
354: case t_QUAD:
355: if (gcmp0((GEN)x[3])) return gtolong((GEN)x[2]); break;
356: }
357: err(typeer,"gtolong");
358: return 0; /* not reached */
359: }
360:
361: /*******************************************************************/
362: /* */
363: /* COMPARE TO ZERO */
364: /* returns 1 whenever the GEN x is 0, and 0 otherwise */
365: /* */
366: /*******************************************************************/
367:
368: int
369: gcmp0(GEN x)
370: {
371: switch(typ(x))
372: {
373: case t_INT: case t_REAL: case t_POL: case t_SER:
374: return !signe(x);
375:
376: case t_INTMOD: case t_POLMOD:
377: return gcmp0((GEN)x[2]);
378:
379: case t_FRAC: case t_FRACN:
380: return !signe(x[1]);
381:
382: case t_COMPLEX:
383: /* is 0 iff norm(x) would be 0 (can happen with Re(x) and Im(x) != 0
384: * only if Re(x) and Im(x) are of type t_REAL). See mp.c:addrr().
385: */
386: if (gcmp0((GEN)x[1]))
387: {
388: if (gcmp0((GEN)x[2])) return 1;
389: if (typ(x[1])!=t_REAL || typ(x[2])!=t_REAL) return 0;
390: return (expo(x[1])>expo(x[2]));
391: }
392: if (gcmp0((GEN)x[2]))
393: {
394: if (typ(x[1])!=t_REAL || typ(x[2])!=t_REAL) return 0;
395: return (expo(x[2])>expo(x[1]));
396: }
397: return 0;
398:
399: case t_PADIC:
400: return !signe(x[4]);
401:
402: case t_QUAD:
403: return gcmp0((GEN)x[2]) && gcmp0((GEN)x[3]);
404:
405: case t_RFRAC: case t_RFRACN:
406: return gcmp0((GEN)x[1]);
407:
408: case t_VEC: case t_COL: case t_MAT:
409: {
410: long i;
411: for (i=lg(x)-1; i; i--)
412: if (!gcmp0((GEN)x[i])) return 0;
413: return 1;
414: }
415: }
416: return 0;
417: }
418:
419: /*******************************************************************/
420: /* */
421: /* COMPARE TO 1 and -1 */
422: /* returns 1 whenever the GEN x is 1 (resp. -1), 0 otherwise */
423: /* */
424: /*******************************************************************/
425:
426: int
427: gcmp1(GEN x)
428: {
429: switch(typ(x))
430: {
431: case t_INT:
432: return is_pm1(x) && signe(x)==1;
433:
434: case t_REAL:
435: if (signe(x) > 0 && expo(x)==0 && x[2]==HIGHBIT)
436: {
437: long i,lx = lg(x);
438: for (i=3; i<lx; i++)
439: if (x[i]) return 0;
440: return 1;
441: }
442: return 0;
443:
444: case t_INTMOD: case t_POLMOD:
445: return gcmp1((GEN)x[2]);
446:
447: case t_FRAC:
448: return gcmp1((GEN)x[1]) && gcmp1((GEN)x[2]);
449:
450: case t_FRACN:
451: return egalii((GEN)x[1],(GEN)x[2]);
452:
453: case t_COMPLEX:
454: return gcmp1((GEN)x[1]) && gcmp0((GEN)x[2]);
455:
456: case t_PADIC:
457: return !valp(x) && gcmp1((GEN)x[4]);
458:
459: case t_QUAD:
460: return gcmp1((GEN)x[2]) && gcmp0((GEN)x[3]);
461:
462: case t_POL:
463: return lgef(x)==3 && gcmp1((GEN)x[2]);
464: }
465: return 0;
466: }
467:
468: int
469: gcmp_1(GEN x)
470: {
471: long l,y;
472: GEN p1;
473:
474: switch(typ(x))
475: {
476: case t_INT:
477: return is_pm1(x) && signe(x)== -1;
478:
479: case t_REAL:
480: if (signe(x) < 0 && expo(x)==0 && x[2]==HIGHBIT)
481: {
482: long i,lx = lg(x);
483: for (i=3; i<lx; i++)
484: if (x[i]) return 0;
485: return 1;
486: }
487: return 0;
488:
489: case t_INTMOD:
490: l=avma; y=egalii(addsi(1,(GEN)x[2]), (GEN)x[1]); avma=l; return y;
491:
492: case t_FRAC: case t_FRACN:
493: l = signe(x[1]);
494: return l && l == -signe(x[2]) && !absi_cmp((GEN)x[1],(GEN)x[2]);
495:
496: case t_COMPLEX:
497: return gcmp_1((GEN)x[1]) && gcmp0((GEN)x[2]);
498:
499: case t_QUAD:
500: return gcmp_1((GEN)x[2]) && gcmp0((GEN)x[3]);
501:
502: case t_PADIC:
503: l=avma; y=gegal(addsi(1,(GEN)x[4]), (GEN)x[3]); avma=l; return y;
504:
505: case t_POLMOD:
506: l=avma; p1=gadd(gun,(GEN)x[2]);
507: y = signe(p1) && !gegal(p1,(GEN)x[1]); avma=l; return !y;
508:
509: case t_POL:
510: return lgef(x)==3 && gcmp_1((GEN)x[2]);
511: }
512: return 0;
513: }
514:
515: /*******************************************************************/
516: /* */
517: /* SIGNED COMPARISON */
518: /* returns the sign of x - y when it makes sense. 0 otherwise */
519: /* */
520: /*******************************************************************/
521:
522: int
523: gcmp(GEN x, GEN y)
524: {
525: long tx,ty,f,av;
526:
527: tx=typ(x); ty=typ(y);
528: if (is_intreal_t(tx))
529: { if (is_intreal_t(ty)) return mpcmp(x,y); }
530: else
531: {
532: if (tx==t_STR)
533: {
534: if (ty != t_STR) return 1;
535: return strcmp(GSTR(x),GSTR(y));
536: }
537: if (!is_frac_t(tx)) err(typeer,"comparison");
538: }
539: if (ty == t_STR) return -1;
540: if (!is_intreal_t(ty) && !is_frac_t(ty)) err(typeer,"comparison");
541: av=avma; y=gneg_i(y); f=gsigne(gadd(x,y)); avma=av; return f;
542: }
543:
544: int
545: lexcmp(GEN x, GEN y)
546: {
547: const long tx=typ(x), ty=typ(y);
548: long lx,ly,l,fl,i;
549:
550: ly=lg(y);
551: if (!is_matvec_t(tx))
552: {
553: if (!is_matvec_t(ty)) return gcmp(x,y);
554: if (ly==1) return 1;
555: fl = lexcmp(x,(GEN)y[1]);
556: if (fl) return fl;
557: return (ly>2)? -1:0;
558: }
559:
560: lx=lg(x);
561: if (!is_matvec_t(ty))
562: {
563: if (lx==1) return -1;
564: fl = lexcmp(y,(GEN)x[1]);
565: if (fl) return -fl;
566: return (lx>2)? 1:0;
567: }
568:
569: /* x and y are matvec_t */
570: if (ly==1) return (lx==1)?0:1;
571: if (lx==1) return -1;
572: if (ty==t_MAT)
573: {
574: if (tx != t_MAT)
575: {
576: fl = lexcmp(x,(GEN)y[1]);
577: if (fl) return fl;
578: return (ly>2)?-1:0;
579: }
580: }
581: else if (tx==t_MAT)
582: {
583: fl = lexcmp(y,(GEN)x[1]);
584: if (fl) return -fl;
585: return (ly>2)?1:0;
586: }
587:
588: /* tx = ty = t_MAT, or x and y are both vect_t */
589: l=min(lx,ly);
590: for (i=1; i<l; i++)
591: {
592: fl = lexcmp((GEN)x[i],(GEN)y[i]);
593: if (fl) return fl;
594: }
595: return (ly != lx)? -1 : 0;
596: }
597:
598: /*****************************************************************/
599: /* */
600: /* EQUALITY */
601: /* returns 1 if x == y, 0 otherwise */
602: /* */
603: /*****************************************************************/
604:
605: static int
606: polegal(GEN x, GEN y)
607: {
608: long i, lx = lgef(x);
609:
610: if (x[1] != y[1] && (lgef(y) != lx || lx > 3)) return 0;
611: for (i = 2; i < lx; i++)
612: if (!gegal((GEN)x[i],(GEN)y[i])) return 0;
613: return 1;
614: }
615:
616: #define MASK(x) (((ulong)(x)) & (TYPBITS | LGBITS))
617: static int
618: vecegal(GEN x, GEN y)
619: {
620: long i, tx = typ(x);
621:
622: if (!is_matvec_t(tx)) return gegal(x,y);
623: if (MASK(x[0]) != MASK(y[0])) return 0;
624:
625: i = lg(x)-1;
626: if (tx != t_MAT)
627: {
628: for ( ; i; i--)
629: if (! gegal((GEN)x[i],(GEN)y[i]) ) return 0;
630: return 1;
631: }
632: for ( ; i; i--)
633: if (! vecegal((GEN)x[i],(GEN)y[i]) ) return 0;
634: return 1;
635: }
636: #undef MASK
637:
638: #define MASK(x) (((ulong)(x)) & (LGEFINTBITS | SIGNBITS))
639: int
640: egalii(GEN x, GEN y)
641: {
642: long i;
643: if (MASK(x[1]) != MASK(y[1])) return 0;
644: i = lgefint(x)-1; while (i>1 && x[i]==y[i]) i--;
645: return i==1;
646: }
647: #undef MASK
648:
649: int
650: gegal(GEN x, GEN y)
651: {
652: long av,i,tx,ty;
653:
654: if (x == y) return 1;
655: tx = typ(x); ty = typ(y);
656: if (tx!=ty)
657: { if (tx == t_STR || ty == t_STR) return 0; }
658: else
659: switch(tx)
660: {
661: case t_INT:
662: return egalii(x,y);
663:
664: case t_POL:
665: return polegal(x,y);
666:
667: case t_COMPLEX:
668: return gegal((GEN)x[1],(GEN)y[1]) && gegal((GEN)x[2],(GEN)y[2]);
669:
670: case t_INTMOD: case t_POLMOD:
671: return gegal((GEN)x[2],(GEN)y[2])
672: && (x[1]==y[1] || gegal((GEN)x[1],(GEN)y[1]));
673:
674: case t_QFR:
675: if (!gegal((GEN)x[4],(GEN)y[4])) return 0; /* fall through */
676: case t_QUAD: case t_QFI:
677: return gegal((GEN)x[1],(GEN)y[1])
678: && gegal((GEN)x[2],(GEN)y[2])
679: && gegal((GEN)x[3],(GEN)y[3]);
680:
681: case t_FRAC:
682: return gegal((GEN)x[1], (GEN)y[1])
683: && gegal((GEN)x[2], (GEN)y[2]);
684:
685: case t_FRACN: case t_RFRAC: case t_RFRACN:
686: av=avma; i=gegal(gmul((GEN)x[1],(GEN)y[2]),gmul((GEN)x[2],(GEN)y[1]));
687: avma=av; return i;
688:
689: case t_STR:
690: return !strcmp(GSTR(x),GSTR(y));
691:
692: case t_VEC: case t_COL: case t_MAT:
693: return vecegal(x,y);
694: }
695: av=avma; y=gneg_i(y); i=gcmp0(gadd(x,y));
696: avma=av; return i;
697: }
698:
699: /*******************************************************************/
700: /* */
701: /* VALUATION */
702: /* p is either an int or a polynomial. */
703: /* returns the largest exponent of p dividing x when this makes */
704: /* sense : error for types real, integermod and polymod if p does */
705: /* not divide the modulus, q-adic if q!=p. */
706: /* */
707: /*******************************************************************/
708:
709: static long
710: minval(GEN x, GEN p, long first, long lx)
711: {
712: long i,k, val = VERYBIGINT;
713: for (i=first; i<lx; i++)
714: {
715: k = ggval((GEN)x[i],p);
716: if (k < val) val = k;
717: }
718: return val;
719: }
720:
721: long
722: ggval(GEN x, GEN p)
723: {
724: long tx=typ(x), tp=typ(p), av, limit,vx,v,i,val;
725: GEN p1,p2;
726:
727: if (isexactzero(x)) return VERYBIGINT;
728: if (is_const_t(tx) && tp==t_POL) return 0;
729: switch(tx)
730: {
731: case t_INT:
732: if (tp != t_INT) break;
733: av=avma; val = pvaluation(x,p, &p1);
734: avma=av; return val;
735:
736: case t_INTMOD:
737: av=avma;
738: p1=cgeti(lgefint(x[1]));
739: p2=cgeti(lgefint(x[2]));
740: if (tp!=t_INT || !mpdivis((GEN)x[1],p,p1)) break;
741: if (!mpdivis((GEN)x[2],p,p2)) { avma=av; return 0; }
742: val=1; while (mpdivis(p1,p,p1) && mpdivis(p2,p,p2)) val++;
743: avma=av; return val;
744:
745: case t_PADIC:
746: if (tp!=t_INT || !gegal(p,(GEN)x[2])) break;
747: return valp(x);
748:
749: case t_POLMOD:
750: if (tp==t_INT) return ggval((GEN)x[2],p);
751: av = avma;
752: if (tp!=t_POL || !poldivis((GEN)x[1],p,&p1)) break;
753: if (!poldivis((GEN)x[2],p,&p2)) { avma=av; return 0; }
754: val=1; while (poldivis(p1,p,&p1)&&poldivis(p2,p,&p2)) val++;
755: avma = av; return val;
756:
757: case t_POL:
758: if (tp==t_POL)
759: {
760: v=varn(p); vx=varn(x);
761: if (vx == v)
762: {
763: if ((p>=(GEN)polx && p <= (GEN)(polx+MAXVARN)) || ismonome(p))
764: {
765: i=2; while (isexactzero((GEN)x[i])) i++;
766: return i-2;
767: }
768: av = avma; limit=stack_lim(av,1);
769: for (val=0; ; val++)
770: {
771: if (!poldivis(x,p,&x)) break;
772: if (low_stack(limit, stack_lim(av,1)))
773: {
774: if(DEBUGMEM>1) err(warnmem,"ggval");
775: x = gerepileupto(av, gcopy(x));
776: }
777: }
778: avma = av; return val;
779: }
780: if (vx > v) return 0;
781: }
782: else
783: {
784: if (tp!=t_INT) break;
785: i=2; while (isexactzero((GEN)x[i])) i++;
786: }
787: return minval(x,p,i,lgef(x));
788:
789: case t_SER:
790: if (tp!=t_POL && tp!=t_SER && tp!=t_INT) break;
791: v=gvar(p); vx=varn(x);
792: if (vx==v) return (long)(valp(x)/ggval(p,polx[v]));
793: if (vx>v) return 0;
794: return minval(x,p,2,lg(x));
795:
796: case t_FRAC: case t_FRACN: case t_RFRAC: case t_RFRACN:
797: return ggval((GEN)x[1],p) - ggval((GEN)x[2],p);
798:
799: case t_COMPLEX: case t_QUAD: case t_VEC: case t_COL: case t_MAT:
800: return minval(x,p,1,lg(x));
801: }
802: err(talker,"forbidden or conflicting type in gval");
803: return 0; /* not reached */
804: }
805:
806: /* x is non-zero */
807: long
808: svaluation(ulong x, ulong p, long *py)
809: {
810: ulong v = 0;
811: for(;;)
812: {
813: if (x%p) { *py = x; return v; }
814: v++; x/=p;
815: }
816: }
817:
818: /* x is a non-zero integer */
819: long
820: pvaluation(GEN x, GEN p, GEN *py)
821: {
822: long av,v;
823: GEN p1,p2;
824:
825: if (!is_bigint(x))
826: {
827: long y;
828: if (!is_bigint(p))
829: {
830: v = svaluation(x[2],p[2], &y);
831: if (signe(x) < 0) y = -y;
832: *py = stoi(y);
833: }
834: else
835: {
836: v = 0;
837: *py = icopy(x);
838: }
839: return v;
840: }
841: av = avma; v = 0; (void)new_chunk(lgefint(x));
842: for(;;)
843: {
844: p1 = dvmdii(x,p,&p2);
845: if (p2 != gzero) { avma=av; *py = icopy(x); return v; }
846: v++; x = p1;
847: }
848: }
849:
850: /*******************************************************************/
851: /* */
852: /* NEGATION: Create -x */
853: /* */
854: /*******************************************************************/
855:
856: GEN
857: gneg(GEN x)
858: {
859: long tx=typ(x),lx,i;
860: GEN y;
861:
862: if (gcmp0(x)) return gcopy(x);
863: switch(tx)
864: {
865: case t_INT: case t_REAL:
866: return mpneg(x);
867:
868: case t_INTMOD: y=cgetg(3,t_INTMOD);
869: icopyifstack(x[1],y[1]);
870: y[2] = lsubii((GEN)y[1],(GEN)x[2]);
871: break;
872:
873: case t_POLMOD: y=cgetg(3,t_POLMOD);
874: copyifstack(x[1],y[1]);
875: y[2]=lneg((GEN)x[2]); break;
876:
877: case t_FRAC: case t_FRACN: case t_RFRAC: case t_RFRACN:
878: y=cgetg(3,tx);
879: y[1]=lneg((GEN)x[1]);
880: y[2]=lcopy((GEN)x[2]); break;
881:
882: case t_PADIC:
883: y=cgetp2(x,valp(x));
884: y[4]=lsubii((GEN)x[3],(GEN)x[4]);
885: break;
886:
887: case t_QUAD:
888: y=cgetg(4,t_QUAD); copyifstack(x[1],y[1]);
889: y[2]=lneg((GEN)x[2]);
890: y[3]=lneg((GEN)x[3]); break;
891:
892: case t_COMPLEX: case t_VEC: case t_COL: case t_MAT:
893: lx=lg(x); y=cgetg(lx,tx);
894: for (i=1; i<lx; i++) y[i]=lneg((GEN)x[i]);
895: break;
896:
897: case t_POL:
898: lx=lgef(x); y=cgetg(lx,tx); y[1]=x[1];
899: for (i=2; i<lx; i++) y[i]=lneg((GEN)x[i]);
900: break;
901:
902: case t_SER:
903: lx=lg(x); y=cgetg(lx,tx); y[1]=x[1];
904: for (i=2; i<lx; i++) y[i]=lneg((GEN)x[i]);
905: break;
906:
907: default:
908: err(typeer,"negation");
909: }
910: return y;
911: }
912:
913: GEN
914: gneg_i(GEN x)
915: {
916: long tx=typ(x),lx,i;
917: GEN y;
918:
919: if (gcmp0(x)) return x;
920: switch(tx)
921: {
922: case t_INT: case t_REAL:
923: return mpneg(x);
924:
925: case t_INTMOD: y=cgetg(3,t_INTMOD); y[1]=x[1];
926: y[2] = lsubii((GEN)y[1],(GEN)x[2]);
927: break;
928:
929: case t_FRAC: case t_FRACN: case t_RFRAC: case t_RFRACN:
930: y=cgetg(3,tx); y[2]=x[2];
931: y[1]=(long)gneg_i((GEN)x[1]); break;
932:
933: case t_PADIC: y = cgetg(5,t_PADIC); y[2]=x[2]; y[3]=x[3];
934: y[1] = evalprecp(precp(x)) | evalvalp(valp(x));
935: y[4]=lsubii((GEN)x[3],(GEN)x[4]); break;
936:
937: case t_POLMOD: y=cgetg(3,t_POLMOD); y[1]=x[1];
938: y[2]=(long)gneg_i((GEN)x[2]); break;
939:
940: case t_QUAD: y=cgetg(4,t_QUAD); y[1]=x[1];
941: y[2]=(long)gneg_i((GEN)x[2]);
942: y[3]=(long)gneg_i((GEN)x[3]); break;
943:
944: case t_COMPLEX: case t_VEC: case t_COL: case t_MAT:
945: lx=lg(x); y=cgetg(lx,tx);
946: for (i=1; i<lx; i++) y[i]=(long)gneg_i((GEN)x[i]);
947: break;
948:
949: case t_POL: lx=lgef(x); y=cgetg(lx,tx); y[1]=x[1];
950: for (i=2; i<lx; i++) y[i]=(long)gneg_i((GEN)x[i]);
951: break;
952:
953: case t_SER: lx=lg(x); y=cgetg(lx,tx); y[1]=x[1];
954: for (i=2; i<lx; i++) y[i]=(long)gneg_i((GEN)x[i]);
955: break;
956:
957: default:
958: err(typeer,"negation");
959: }
960: return y;
961: }
962:
963: /******************************************************************/
964: /* */
965: /* ABSOLUTE VALUE */
966: /* Create abs(x) if x is integer, real, fraction or complex. */
967: /* Error otherwise. */
968: /* */
969: /******************************************************************/
970:
971: GEN
972: gabs(GEN x, long prec)
973: {
974: long tx=typ(x),lx,i,l,tetpil;
975: GEN y,p1;
976:
977: switch(tx)
978: {
979: case t_INT: case t_REAL:
980: return mpabs(x);
981:
982: case t_FRAC: case t_FRACN: y=cgetg(lg(x),tx);
983: y[1]=labsi((GEN)x[1]);
984: y[2]=labsi((GEN)x[2]); return y;
985:
986: case t_COMPLEX:
987: l=avma; p1=gnorm(x); tetpil=avma;
988: return gerepile(l,tetpil,gsqrt(p1,prec));
989:
990: case t_QUAD:
991: l=avma; p1=gmul(x, realun(prec)); tetpil=avma;
992: return gerepile(l,tetpil,gabs(p1,prec));
993:
994: case t_POL:
995: lx=lgef(x); if (lx<=2) return gcopy(x);
996: p1=(GEN)x[lx-1];
997: switch(typ(p1))
998: {
999: case t_INT: case t_REAL: case t_FRAC: case t_FRACN:
1000: if (gsigne(p1)<0) return gneg(x);
1001: }
1002: return gcopy(x);
1003:
1004: case t_VEC: case t_COL: case t_MAT:
1005: lx=lg(x); y=cgetg(lx,tx);
1006: for (i=1; i<lx; i++)
1007: y[i]=(long)gabs((GEN)x[i],prec);
1008: return y;
1009: }
1010: err(typeer,"gabs");
1011: return NULL; /* not reached */
1012: }
1013:
1014: GEN
1015: gmax(GEN x, GEN y)
1016: {
1017: if (gcmp(x,y)>0) y=x; return gcopy(y);
1018: }
1019:
1020: GEN
1021: gmin(GEN x, GEN y)
1022: {
1023: if (gcmp(x,y)<0) y=x; return gcopy(y);
1024: }
1025:
1026: GEN
1027: vecmax(GEN x)
1028: {
1029: long tx=typ(x),lx,lx2,i,j;
1030: GEN *p1,s;
1031:
1032: if (!is_matvec_t(tx)) return gcopy(x);
1033: lx=lg(x); if (lx==1) return stoi(-VERYBIGINT);
1034: if (tx!=t_MAT)
1035: {
1036: s=(GEN)x[1];
1037: for (i=2; i<lx; i++)
1038: if (gcmp((GEN)x[i],s) > 0) s=(GEN)x[i];
1039: }
1040: else
1041: {
1042: lx2 = lg(x[1]);
1043: if (lx2==1) return stoi(-VERYBIGINT);
1044: s=gmael(x,1,1); i=2;
1045: for (j=1; j<lx; j++)
1046: {
1047: p1 = (GEN *) x[j];
1048: for (; i<lx2; i++)
1049: if (gcmp(p1[i],s) > 0) s=p1[i];
1050: i=1;
1051: }
1052: }
1053: return gcopy(s);
1054: }
1055:
1056: GEN
1057: vecmin(GEN x)
1058: {
1059: long tx=typ(x),lx,lx2,i,j;
1060: GEN *p1,s;
1061:
1062: if (!is_matvec_t(tx)) return gcopy(x);
1063: lx=lg(x); if (lx==1) return stoi(VERYBIGINT);
1064: if (tx!=t_MAT)
1065: {
1066: s=(GEN)x[1];
1067: for (i=2; i<lx; i++)
1068: if (gcmp((GEN)x[i],s) < 0) s=(GEN)x[i];
1069: }
1070: else
1071: {
1072: lx2 = lg(x[1]);
1073: if (lx2==1) return stoi(VERYBIGINT);
1074: s=gmael(x,1,1); i=2;
1075: for (j=1; j<lx; j++)
1076: {
1077: p1 = (GEN *) x[j];
1078: for (; i<lx2; i++)
1079: if (gcmp(p1[i],s) < 0) s=p1[i];
1080: i=1;
1081: }
1082: }
1083: return gcopy(s);
1084: }
1085:
1086: /*******************************************************************/
1087: /* */
1088: /* AFFECT long --> GEN */
1089: /* affect long s to GEN x. Useful for initialization. */
1090: /* */
1091: /*******************************************************************/
1092:
1093: static void
1094: padicaff0(GEN x)
1095: {
1096: if (signe(x[4]))
1097: {
1098: setvalp(x, valp(x)|precp(x));
1099: affsi(0,(GEN)x[4]);
1100: }
1101: }
1102:
1103: void
1104: gaffsg(long s, GEN x)
1105: {
1106: long l,i,v;
1107: GEN p;
1108:
1109: switch(typ(x))
1110: {
1111: case t_INT:
1112: affsi(s,x); break;
1113:
1114: case t_REAL:
1115: affsr(s,x); break;
1116:
1117: case t_INTMOD:
1118: modsiz(s,(GEN)x[1],(GEN)x[2]); break;
1119:
1120: case t_FRAC: case t_FRACN:
1121: affsi(s,(GEN)x[1]); affsi(1,(GEN)x[2]); break;
1122:
1123: case t_COMPLEX:
1124: gaffsg(s,(GEN)x[1]); gaffsg(0,(GEN)x[2]); break;
1125:
1126: case t_PADIC:
1127: if (!s) { padicaff0(x); break; }
1128: p = (GEN)x[2];
1129: v = (is_bigint(p))? 0: svaluation(s,p[2], &i);
1130: setvalp(x,v); modsiz(i,(GEN)x[3],(GEN)x[4]);
1131: break;
1132:
1133: case t_QUAD:
1134: gaffsg(s,(GEN)x[2]); gaffsg(0,(GEN)x[3]); break;
1135:
1136: case t_POLMOD:
1137: gaffsg(s,(GEN)x[2]); break;
1138:
1139: case t_POL:
1140: v=varn(x);
1141: if (!s) x[1]=evallgef(2) | evalvarn(v);
1142: else
1143: {
1144: x[1]=evalsigne(1) | evallgef(3) | evalvarn(v);
1145: gaffsg(s,(GEN)x[2]);
1146: }
1147: break;
1148:
1149: case t_SER:
1150: v=varn(x); gaffsg(s,(GEN)x[2]); l=lg(x);
1151: if (!s)
1152: x[1] = evalvalp(l-2) | evalvarn(v);
1153: else
1154: x[1] = evalsigne(1) | evalvalp(0) | evalvarn(v);
1155: for (i=3; i<l; i++) gaffsg(0,(GEN)x[i]);
1156: break;
1157:
1158: case t_RFRAC: case t_RFRACN:
1159: gaffsg(s,(GEN)x[1]); gaffsg(1,(GEN)x[2]); break;
1160:
1161: case t_VEC: case t_COL: case t_MAT:
1162: if (lg(x)!=2) err(assigneri,t_INT,typ(x));
1163: gaffsg(s,(GEN)x[1]); break;
1164:
1165: default: err(typeer,"gaffsg");
1166: }
1167: }
1168:
1169: /*******************************************************************/
1170: /* */
1171: /* GENERIC AFFECTATION */
1172: /* Affect the content of x to y, whenever possible */
1173: /* */
1174: /*******************************************************************/
1175:
1176: GEN
1177: realzero(long prec)
1178: {
1179: GEN x=cgetr(3);
1180: x[1]=evalexpo(-bit_accuracy(prec));
1181: x[2]=0; return x;
1182: }
1183: GEN
1184: realun(long prec)
1185: {
1186: GEN x=cgetr(prec); affsr(1,x);
1187: return x;
1188: }
1189:
1190: void
1191: gaffect(GEN x, GEN y)
1192: {
1193: long i,j,k,l,v,vy,lx,ly,tx,ty,av;
1194: GEN p1,num,den;
1195:
1196:
1197: tx=typ(x); ty=typ(y); lx=lg(x); ly=lg(y);
1198: if (is_scalar_t(tx))
1199: {
1200: if (is_scalar_t(ty))
1201: {
1202: switch(tx)
1203: {
1204: case t_INT:
1205: switch(ty)
1206: {
1207: case t_INT:
1208: /* y = gzero, gnil, gun or gdeux */
1209: if (is_universal_constant(y))
1210: {
1211: if (y==gzero) err(overwriter,"integer: 0 (gzero)");
1212: if (y==gun) err(overwriter,"integer: 1 (gun)");
1213: if (y==gdeux) err(overwriter,"integer: 2 (gdeux)");
1214: err(overwriter,"integer: nil (gnil)");
1215: }
1216: affii(x,y); break;
1217:
1218: case t_REAL:
1219: if (y == gpi || y==geuler) err(overwriter,"real");
1220: affir(x,y); break;
1221:
1222: case t_INTMOD:
1223: modiiz(x,(GEN)y[1],(GEN)y[2]); break;
1224:
1225: case t_FRAC:
1226: if (y == ghalf) err(overwriter,"fraction: 1/2 (ghalf)");
1227: case t_FRACN:
1228: affii(x,(GEN)y[1]); affsi(1,(GEN)y[2]); break;
1229:
1230: case t_COMPLEX:
1231: if (y == gi) err(overwriter,"complex number: I (gi)");
1232: gaffect(x,(GEN)y[1]); gaffsg(0,(GEN)y[2]); break;
1233:
1234: case t_PADIC:
1235: if (!signe(x)) { padicaff0(y); break; }
1236: l=avma;
1237: setvalp(y, pvaluation(x,(GEN)y[2],&p1));
1238: modiiz(p1,(GEN)y[3],(GEN)y[4]);
1239: avma=l; break;
1240:
1241: case t_QUAD:
1242: gaffect(x,(GEN)y[2]); gaffsg(0,(GEN)y[3]); break;
1243:
1244: case t_POLMOD:
1245: gaffect(x,(GEN)y[2]); break;
1246:
1247: default: err(typeer,"gaffect");
1248: }
1249: break;
1250:
1251: case t_REAL:
1252: switch(ty)
1253: {
1254: case t_REAL:
1255: affrr(x,y); break;
1256:
1257: case t_COMPLEX:
1258: gaffect(x,(GEN)y[1]); gaffsg(0,(GEN)y[2]); break;
1259:
1260: case t_POLMOD:
1261: gaffect(x,(GEN)y[2]); break;
1262:
1263: case t_INT: case t_INTMOD: case t_FRAC:
1264: case t_FRACN: case t_PADIC: case t_QUAD:
1265: err(assignerf,tx,ty);
1266:
1267: default: err(typeer,"gaffect");
1268: }
1269: break;
1270:
1271: case t_INTMOD:
1272: switch(ty)
1273: {
1274: case t_INTMOD:
1275: if (!divise((GEN)x[1],(GEN)y[1]))
1276: err(assigneri,tx,ty);
1277: modiiz((GEN)x[2],(GEN)y[1],(GEN)y[2]); break;
1278:
1279: case t_POLMOD:
1280: gaffect(x,(GEN)y[2]); break;
1281:
1282: case t_INT: case t_REAL: case t_FRAC:
1283: case t_FRACN: case t_COMPLEX: case t_QUAD:
1284: err(assignerf,tx,ty);
1285:
1286: case t_PADIC:
1287: err(assignerf,tx,ty);
1288: default: err(typeer,"gaffect");
1289: }
1290: break;
1291:
1292: case t_FRAC:
1293: switch(ty)
1294: {
1295: case t_INT:
1296: if (! mpdivis((GEN)x[1],(GEN)x[2],y))
1297: err(assigneri,tx,ty);
1298: break;
1299:
1300: case t_REAL:
1301: diviiz((GEN)x[1],(GEN)x[2],y); break;
1302:
1303: case t_INTMOD:
1304: av=avma; p1=mpinvmod((GEN)x[2],(GEN)y[1]);
1305: modiiz(mulii((GEN)x[1],p1),(GEN)y[1],(GEN)y[2]);
1306: avma=av; break;
1307:
1308: case t_FRAC: case t_FRACN:
1309: affii((GEN)x[1],(GEN)y[1]);
1310: affii((GEN)x[2],(GEN)y[2]);
1311: break;
1312:
1313: case t_COMPLEX:
1314: gaffect(x,(GEN)y[1]); gaffsg(0,(GEN)y[2]); break;
1315:
1316: case t_PADIC:
1317: if (!signe(x[1])) { padicaff0(y); break; }
1318: av=avma; num=(GEN)x[1]; den=(GEN)x[2];
1319: v = pvaluation(num, (GEN) y[2], &num);
1320: if (!v) v = -pvaluation(den,(GEN)y[2],&den);
1321: p1=mulii(num,mpinvmod(den,(GEN)y[3]));
1322: modiiz(p1,(GEN)y[3],(GEN)y[4]); avma=av;
1323: setvalp(y,v); break;
1324:
1325: case t_QUAD:
1326: gaffect(x,(GEN)y[2]); gaffsg(0,(GEN)y[3]); break;
1327: case t_POLMOD:
1328: gaffect(x,(GEN)y[2]); break;
1329: default: err(typeer,"gaffect");
1330: }
1331: break;
1332:
1333: case t_FRACN:
1334: switch(ty)
1335: {
1336: case t_INT:
1337: if (! mpdivis((GEN)x[1],(GEN)x[2],y)) err(assigneri,tx,ty);
1338: break;
1339:
1340: case t_REAL:
1341: diviiz((GEN)x[1],(GEN)x[2],y); break;
1342: case t_INTMOD:
1343: av=avma; gaffect(gred(x),y); avma=av; break;
1344: case t_FRAC:
1345: gredz(x,y); break;
1346: case t_FRACN:
1347: affii((GEN)x[1],(GEN)y[1]);
1348: affii((GEN)x[2],(GEN)y[2]); break;
1349: case t_COMPLEX:
1350: gaffect(x,(GEN)y[1]); gaffsg(0,(GEN)y[2]); break;
1351: case t_PADIC:
1352: if (!signe(x[1])) { padicaff0(y); break; }
1353: av=avma; num=(GEN)x[1]; den=(GEN)x[2];
1354: v = pvaluation(num,(GEN)y[2],&num)
1355: - pvaluation(den,(GEN)y[2],&den);
1356: p1=mulii(num,mpinvmod(den,(GEN)y[3]));
1357: modiiz(p1,(GEN)y[3],(GEN)y[4]); avma=av;
1358: setvalp(y,v); break;
1359:
1360: case t_QUAD:
1361: gaffect(x,(GEN)y[2]); gaffsg(0,(GEN)y[3]); break;
1362: case t_POLMOD:
1363: gaffect(x,(GEN)y[2]); break;
1364: default: err(typeer,"gaffect");
1365: }
1366: break;
1367:
1368: case t_COMPLEX:
1369: switch(ty)
1370: {
1371: case t_INT: case t_REAL: case t_INTMOD:
1372: case t_FRAC: case t_FRACN: case t_PADIC: case t_QUAD:
1373: if (!gcmp0((GEN)x[2])) err(assigneri,tx,ty);
1374: gaffect((GEN)x[1],y); break;
1375:
1376: case t_COMPLEX:
1377: gaffect((GEN)x[1],(GEN)y[1]);
1378: gaffect((GEN)x[2],(GEN)y[2]);
1379: break;
1380:
1381: case t_POLMOD:
1382: gaffect(x,(GEN)y[2]); break;
1383:
1384: default: err(typeer,"gaffect");
1385: }
1386: break;
1387:
1388: case t_PADIC:
1389: switch(ty)
1390: {
1391: case t_INTMOD:
1392: if (valp(x)<0) err(assigneri,tx,ty);
1393: av=avma;
1394: v = pvaluation((GEN)y[1],(GEN)x[2],&p1);
1395: k = signe(x[4])? (precp(x)+valp(x)): valp(x)+1;
1396: if (!gcmp1(p1) || v > k) err(assigneri,tx,ty);
1397: modiiz(gtrunc(x),(GEN)y[1],(GEN)y[2]); avma=av; break;
1398:
1399: case t_PADIC:
1400: if (!egalii((GEN)x[2],(GEN)y[2])) err(assigneri,tx,ty);
1401: modiiz((GEN)x[4],(GEN)y[3],(GEN)y[4]);
1402: setvalp(y,valp(x)); break;
1403:
1404: case t_POLMOD:
1405: gaffect(x,(GEN)y[2]); break;
1406:
1407: case t_INT: case t_REAL: case t_FRAC:
1408: case t_FRACN: case t_COMPLEX: case t_QUAD:
1409: err(assignerf,tx,ty);
1410:
1411: default: err(typeer,"gaffect");
1412: }
1413: break;
1414:
1415: case t_QUAD:
1416: switch(ty)
1417: {
1418: case t_INT: case t_INTMOD: case t_FRAC:
1419: case t_FRACN: case t_PADIC:
1420: if (!gcmp0((GEN)x[3])) err(assigneri,tx,ty);
1421: gaffect((GEN)x[2],y); break;
1422:
1423: case t_REAL:
1424: av=avma; p1=co8(x,ly); gaffect(p1,y); avma=av; break;
1425: case t_COMPLEX:
1426: ly=precision(y);
1427: if (ly)
1428: {
1429: av=avma; p1=co8(x,ly); gaffect(p1,y); avma=av;
1430: }
1431: else
1432: {
1433: if (!gcmp0((GEN)x[3])) err(assigneri,tx,ty);
1434: gaffect((GEN)x[2],y);
1435: }
1436: break;
1437:
1438: case t_QUAD:
1439: if (! gegal((GEN)x[1],(GEN)y[1])) err(assigneri,tx,ty);
1440: affii((GEN)x[2],(GEN)y[2]);
1441: affii((GEN)x[3],(GEN)y[3]);
1442: break;
1443: case t_POLMOD:
1444: gaffect(x,(GEN)y[2]); break;
1445: default: err(typeer,"gaffect");
1446: }
1447: break;
1448:
1449: case t_POLMOD:
1450: if (ty!=t_POLMOD) err(assignerf,tx,ty);
1451: if (! gdivise((GEN)x[1],(GEN)y[1])) err(assigneri,tx,ty);
1452: gmodz((GEN)x[2],(GEN)y[1],(GEN)y[2]); break;
1453: default: err(typeer,"gaffect");
1454: }
1455: return;
1456: }
1457:
1458: /* here y is not scalar */
1459: switch(ty)
1460: {
1461: case t_POL:
1462: v=varn(y);
1463: if (y==polun[v] || y==polx[v])
1464: err(talker,"trying to overwrite a universal polynomial");
1465: gaffect(x,(GEN)y[2]);
1466: for (i=3; i<ly; i++) gaffsg(0,(GEN)y[i]);
1467: if (gcmp0(x)) y[1]=evallgef(2)+evalvarn(v);
1468: else y[1]=evalsigne(1)+evallgef(3)+evalvarn(v);
1469: break;
1470:
1471: case t_SER:
1472: v=varn(y); gaffect(x,(GEN)y[2]);
1473: if (gcmp0(x))
1474: y[1] = evalvalp(ly-2) | evalvarn(v);
1475: else
1476: y[1] = evalsigne(1) | evalvalp(0) | evalvarn(v);
1477: for (i=3; i<ly; i++) gaffsg(0,(GEN)y[i]);
1478: break;
1479:
1480: case t_RFRAC: case t_RFRACN:
1481: gaffect(x,(GEN)y[1]); gaffsg(1,(GEN)y[2]); break;
1482:
1483: case t_VEC: case t_COL: case t_MAT:
1484: err(assignerf,tx,ty);
1485:
1486: default: err(typeer,"gaffect");
1487: }
1488: return;
1489: }
1490:
1491: if (is_const_t(ty)) {
1492: /* initial_value macro is not defined now... */
1493: #define initial_value(ep) ((ep)+1)
1494: if (tx == t_POL) {
1495: entree *var = varentries[ordvar[varn(x)]];
1496: /* Is a bound expression, thus should not loop: */
1497: if (var && var->value != (void*)initial_value(var)) {
1498: gaffect(geval(x),y);
1499: return;
1500: }
1501: } else if (is_rfrac_t(tx)) {
1502: GEN num = (GEN)x[1], den = (GEN)x[2];
1503: entree *varnum = varentries[ordvar[varn(num)]];
1504: entree *varden = varentries[ordvar[varn(den)]];
1505:
1506: /* Are bound expressions, thus should not loop: */
1507: if (varnum && varnum->value != (void*)initial_value(varnum)
1508: && varden && varden->value != (void*)initial_value(varden)) {
1509: gaffect(geval(x),y);
1510: return;
1511: }
1512: }
1513: #undef initial_value
1514: err(assignerf,tx,ty);
1515: }
1516:
1517: switch(tx)
1518: {
1519: case t_POL:
1520: v=varn(x);
1521: switch(ty)
1522: {
1523: case t_POL:
1524: vy=varn(y); if (vy>v) err(assignerf,tx,ty);
1525: if (vy==v)
1526: {
1527: l=lgef(x); if (l>ly) err(assigneri,tx,ty);
1528: y[1]=x[1]; for (i=2; i<l; i++) gaffect((GEN)x[i],(GEN)y[i]);
1529: }
1530: else
1531: {
1532: gaffect(x,(GEN)y[2]);
1533: for (i=3; i<ly; i++) gaffsg(0,(GEN)y[i]);
1534: if (signe(x)) y[1]=evalsigne(1)+evallgef(3)+evalvarn(vy);
1535: else y[1]=evallgef(2)+evalvarn(vy);
1536: }
1537: break;
1538:
1539: case t_SER:
1540: vy=varn(y); if (vy>v) err(assignerf,tx,ty);
1541: if (!signe(x)) { gaffsg(0,y); return; }
1542: if (vy==v)
1543: {
1544: i=gval(x,v); y[1]=evalvarn(v) | evalvalp(i) | evalsigne(1);
1545: k=lgef(x)-i; if (k>ly) k=ly;
1546: for (j=2; j<k; j++) gaffect((GEN)x[i+j],(GEN)y[j]);
1547: for ( ; j<ly; j++) gaffsg(0,(GEN)y[j]);
1548: }
1549: else
1550: {
1551: gaffect(x,(GEN)y[2]);
1552: if (!signe(x))
1553: y[1] = evalvalp(ly-2) | evalvarn(vy);
1554: else
1555: y[1] = evalsigne(1) | evalvalp(0) | evalvarn(vy);
1556: for (i=3; i<ly; i++) gaffsg(0,(GEN)y[i]);
1557: }
1558: break;
1559:
1560: case t_POLMOD:
1561: gmodz(x,(GEN)y[1],(GEN)y[2]); break;
1562:
1563: case t_RFRAC: case t_RFRACN:
1564: gaffect(x,(GEN)y[1]); gaffsg(1,(GEN)y[2]); break;
1565:
1566: case t_VEC: case t_COL: case t_MAT:
1567: err(assignerf,tx,ty);
1568:
1569: default: err(typeer,"gaffect");
1570: }
1571: break;
1572:
1573: case t_SER:
1574: if (ty!=t_SER) err(assignerf,tx,ty);
1575: v=varn(x); vy=varn(y); if (vy>v) err(assignerf,tx,ty);
1576: if (vy==v)
1577: {
1578: y[1]=x[1]; k=lx; if (k>ly) k=ly;
1579: for (i=2; i< k; i++) gaffect((GEN)x[i],(GEN)y[i]);
1580: for ( ; i<ly; i++) gaffsg(0,(GEN)y[i]);
1581: }
1582: else
1583: {
1584: gaffect(x,(GEN)y[2]);
1585: if (!signe(x))
1586: y[1] = evalvalp(ly-2) | evalvarn(vy);
1587: else
1588: y[1] = evalsigne(1) | evalvalp(0) | evalvarn(vy);
1589: for (i=3; i<ly; i++) gaffsg(0,(GEN)y[i]);
1590: }
1591: break;
1592:
1593: case t_RFRAC: case t_RFRACN:
1594: switch(ty)
1595: {
1596: case t_POL: case t_VEC: case t_COL: case t_MAT:
1597: err(assignerf,tx,ty);
1598:
1599: case t_POLMOD:
1600: av=avma; p1=ginvmod((GEN)x[2],(GEN)y[1]);
1601: gmodz(gmul((GEN)x[1],p1),(GEN)y[1],(GEN)y[2]);
1602: avma=av; break;
1603:
1604: case t_SER:
1605: gdivz((GEN)x[1],(GEN)x[2],y); break;
1606:
1607: case t_RFRAC:
1608: if (tx==t_RFRACN) gredz(x,y);
1609: else
1610: {
1611: gaffect((GEN)x[1],(GEN)y[1]);
1612: gaffect((GEN)x[2],(GEN)y[2]);
1613: }
1614: break;
1615:
1616: case t_RFRACN:
1617: gaffect((GEN)x[1],(GEN)y[1]);
1618: gaffect((GEN)x[2],(GEN)y[2]); break;
1619:
1620: default: err(typeer,"gaffect");
1621: }
1622: break;
1623:
1624: case t_QFR: case t_QFI: case t_VEC: case t_COL: case t_MAT:
1625: if (ty != tx || ly != lx) err(assigneri,tx,ty);
1626: for (i=1; i<lx; i++) gaffect((GEN)x[i],(GEN)y[i]);
1627: break;
1628:
1629: default: err(typeer,"gaffect");
1630: }
1631: }
1632:
1633: /*******************************************************************/
1634: /* */
1635: /* CONVERSION QUAD --> REAL, COMPLEX OR P-ADIC */
1636: /* */
1637: /*******************************************************************/
1638:
1639: GEN
1640: co8(GEN x, long prec)
1641: {
1642: long av=avma,tetpil;
1643: GEN p1, p = (GEN) x[1];
1644:
1645: p1 = subii(sqri((GEN)p[3]), shifti((GEN)p[2],2));
1646: if (signe(p1) > 0)
1647: {
1648: p1 = subri(gsqrt(p1,prec), (GEN)p[3]);
1649: setexpo(p1, expo(p1)-1);
1650: }
1651: else
1652: {
1653: p1 = gsub(gsqrt(p1,prec), (GEN)p[3]);
1654: p1[1] = lmul2n((GEN)p1[1],-1);
1655: setexpo(p1[2], expo(p1[2])-1);
1656: }/* p1 = (-b + sqrt(D)) / 2 */
1657: p1 = gmul((GEN)x[3],p1); tetpil=avma;
1658: return gerepile(av,tetpil,gadd((GEN)x[2],p1));
1659: }
1660:
1661: GEN
1662: cvtop(GEN x, GEN p, long l)
1663: {
1664: GEN p1,p2,p3;
1665: long av,tetpil,n;
1666:
1667: if (typ(p)!=t_INT)
1668: err(talker,"not an integer modulus in cvtop or gcvtop");
1669: if (gcmp0(x)) return ggrandocp(p,l);
1670: switch(typ(x))
1671: {
1672: case t_INT:
1673: return gadd(x,ggrandocp(p,ggval(x,p)+l));
1674:
1675: case t_INTMOD:
1676: n=ggval((GEN)x[1],p); if (n>l) n=l;
1677: return gadd((GEN)x[2],ggrandocp(p,n));
1678:
1679: case t_FRAC: case t_FRACN:
1680: n = ggval((GEN)x[1],p) - ggval((GEN)x[2],p);
1681: return gadd(x,ggrandocp(p,n+l));
1682:
1683: case t_COMPLEX:
1684: av=avma; p1=gsqrt(gaddgs(ggrandocp(p,l),-1),0);
1685: p1=gmul(p1,(GEN)x[2]); tetpil=avma;
1686: return gerepile(av,tetpil,gadd(p1,(GEN)x[1]));
1687:
1688: case t_PADIC:
1689: return gprec(x,l);
1690:
1691: case t_QUAD:
1692: av=avma; p1=(GEN)x[1]; p3=gmul2n((GEN)p1[3],-1);
1693: p2=gsub(gsqr(p3),(GEN)p1[2]);
1694: switch(typ(p2))
1695: {
1696: case t_INT:
1697: n=ggval(p2,p); p2=gadd(p2,ggrandocp(p,n+l)); break;
1698:
1699: case t_INTMOD: case t_PADIC:
1700: break;
1701:
1702: case t_FRAC: case t_FRACN:
1703: n = ggval((GEN)p2[1],p) - ggval((GEN)p2[2],p);
1704: p2=gadd(p2,ggrandocp(p,n+l)); break;
1705:
1706: default: err(assigneri,t_QUAD,t_QUAD);
1707: }
1708: p2=gsqrt(p2,0); p1=gmul((GEN)x[3],gsub(p2,p3)); tetpil=avma;
1709: return gerepile(av,tetpil,gadd((GEN)x[2],p1));
1710:
1711: }
1712: err(typeer,"cvtop");
1713: return NULL; /* not reached */
1714: }
1715:
1716: GEN
1717: gcvtop(GEN x, GEN p, long r)
1718: {
1719: long i,lx, tx=typ(x);
1720: GEN y;
1721:
1722: if (is_const_t(tx)) return cvtop(x,p,r);
1723: switch(tx)
1724: {
1725: case t_POL:
1726: lx=lgef(x); y=cgetg(lx,t_POL); y[1]=x[1];
1727: for (i=2; i<lx; i++)
1728: y[i]=(long)gcvtop((GEN)x[i],p,r);
1729: break;
1730:
1731: case t_SER:
1732: lx=lg(x); y=cgetg(lx,t_SER); y[1]=x[1];
1733: for (i=2; i<lx; i++)
1734: y[i]=(long)gcvtop((GEN)x[i],p,r);
1735: break;
1736:
1737: case t_POLMOD: case t_RFRAC: case t_RFRACN:
1738: case t_VEC: case t_COL: case t_MAT:
1739: lx=lg(x); y=cgetg(lx,tx);
1740: for (i=1; i<lx; i++)
1741: y[i]=(long)gcvtop((GEN)x[i],p,r);
1742: break;
1743:
1744: default: err(typeer,"gcvtop");
1745: }
1746: return y;
1747: }
1748:
1749: long
1750: gexpo(GEN x)
1751: {
1752: long tx=typ(x),lx,e,i,y,av;
1753:
1754: switch(tx)
1755: {
1756: case t_INT:
1757: return expi(x);
1758:
1759: case t_FRAC: case t_FRACN:
1760: if (!signe(x[1])) return -HIGHEXPOBIT;
1761: return expi((GEN)x[1]) - expi((GEN)x[2]);
1762:
1763: case t_REAL:
1764: return expo(x);
1765:
1766: case t_COMPLEX:
1767: e = gexpo((GEN)x[1]);
1768: y = gexpo((GEN)x[2]); return max(e,y);
1769:
1770: case t_QUAD:
1771: av=avma; x = co8(x,3); avma=av; return gexpo(x);
1772:
1773: case t_POL: case t_SER: case t_VEC: case t_COL: case t_MAT:
1774: lx=(tx==t_POL)? lgef(x): lg(x);
1775: y = -HIGHEXPOBIT;
1776: for (i=lontyp[tx]; i<lx; i++) { e=gexpo((GEN)x[i]); if (e>y) y=e; }
1777: return y;
1778: }
1779: err(typeer,"gexpo");
1780: return 0; /* not reached */
1781: }
1782:
1783: long
1784: gsize(GEN x)
1785: {
1786: return gcmp0(x)? 0: (long) ((gexpo(x)+1) * L2SL10) + 1;
1787: }
1788:
1789: /* Normalize series x in place.
1790: * Assumption: x,x[2],...,x[lg(x)-1] have been created in that order.
1791: * All intermediate objects will be destroyed.
1792: */
1793: GEN
1794: normalize(GEN x)
1795: {
1796: long i,j, lx = lg(x);
1797:
1798: if (typ(x)!=t_SER) err(typeer,"normalize");
1799: if (lx==2) { setsigne(x,0); avma = (long) x; return x; }
1800: if (! isexactzero((GEN)x[2])) { setsigne(x,1); return x; }
1801:
1802: for (i=3; i<lx; i++)
1803: if (! isexactzero((GEN)x[i]))
1804: {
1805: long tetpil = avma;
1806: GEN p1 = cgetg(lx-i+2,t_SER);
1807: p1[1] = evalsigne(1) | evalvalp(valp(x)+i-2) | evalvarn(varn(x));
1808: j=i; i=2; while (j<lx) p1[i++] = lcopy((GEN)x[j++]);
1809: return gerepile((long) (x+lx),tetpil,p1);
1810: }
1811: avma = (long) (x+lx); return zeroser(varn(x),lx-2);
1812: }
1813:
1814: GEN
1815: normalizepol_i(GEN x, long lx)
1816: {
1817: long i;
1818: for (i=lx-1; i>1; i--)
1819: if (! isexactzero((GEN)x[i])) break;
1820: setlgef(x,i+1);
1821: for (; i>1; i--)
1822: if (! gcmp0((GEN)x[i]) ) { setsigne(x,1); return x; }
1823: setsigne(x,0); return x;
1824: }
1825:
1826: /* Normalize polynomial x in place. See preceding comment */
1827: GEN
1828: normalizepol(GEN x)
1829: {
1830: if (typ(x)!=t_POL) err(typeer,"normalizepol");
1831: return normalizepol_i(x, lgef(x));
1832: }
1833:
1834: int
1835: gsigne(GEN x)
1836: {
1837: switch(typ(x))
1838: {
1839: case t_INT: case t_REAL:
1840: return signe(x);
1841:
1842: case t_FRAC: case t_FRACN:
1843: return (signe(x[2])>0) ? signe(x[1]) : -signe(x[1]);
1844: }
1845: err(typeer,"gsigne");
1846: return 0; /* not reached */
1847: }
1848:
1849: int ff_poltype(GEN *x, GEN *p, GEN *pol);
1850:
1851: GEN
1852: gsqr(GEN x)
1853: {
1854: long tx=typ(x),lx,i,j,k,l,av,tetpil;
1855: GEN z,p1,p2,p3,p4;
1856:
1857: if (is_scalar_t(tx))
1858: switch(tx)
1859: {
1860: case t_INT:
1861: return sqri(x);
1862:
1863: case t_REAL:
1864: return mulrr(x,x);
1865:
1866: case t_INTMOD: z=cgetg(3,t_INTMOD); p2=(GEN)x[1];
1867: (void)new_chunk(lgefint(p2)<<2);
1868: p1=sqri((GEN)x[2]); avma=(long)z;
1869: z[2]=lmodii(p1,p2); icopyifstack(p2,z[1]); return z;
1870:
1871: case t_FRAC: case t_FRACN:
1872: z=cgetg(3,tx);
1873: z[1]=lsqri((GEN)x[1]);
1874: z[2]=lsqri((GEN)x[2]);
1875: return z; /* reduction is useless ! */
1876:
1877: case t_COMPLEX:
1878: z=cgetg(lg(x),tx); l=avma;
1879: p1=gadd((GEN)x[1],(GEN)x[2]);
1880: p2=gadd((GEN)x[1],gneg_i((GEN)x[2]));
1881: p3=gmul((GEN)x[1],(GEN)x[2]);
1882: tetpil=avma;
1883: z[1]=lmul(p1,p2); z[2]=lshift(p3,1);
1884: gerepilemanyvec(l,tetpil,z+1,2);
1885: return z;
1886:
1887: case t_PADIC:
1888: z=cgetg(5,t_PADIC);
1889: z[2] = lcopy((GEN)x[2]);
1890: if (!cmpsi(2,(GEN)x[2]))
1891: {
1892: i=precp(x)+1; av=avma;
1893: p1=addii((GEN)x[3],shifti((GEN)x[4],1));
1894: if (!gcmp0(p1)) { j=vali(p1); if (j<i) i=j; }
1895: avma=av;
1896: z[3]=lshifti((GEN)x[3],i);
1897: z[1]=evalprecp(precp(x)+i) | evalvalp(2*valp(x));
1898: }
1899: else
1900: {
1901: z[3]=licopy((GEN)x[3]);
1902: z[1]=evalprecp(precp(x)) | evalvalp(2*valp(x));
1903: }
1904: z[4]=lgeti(lg(z[3])); av=avma;
1905: modiiz(sqri((GEN)x[4]),(GEN)z[3],(GEN)z[4]);
1906: avma=av; return z;
1907:
1908: case t_QUAD:
1909: p1=(GEN)x[1]; z=cgetg(lg(x),tx); l=avma;
1910: p2=gsqr((GEN)x[2]); p3=gsqr((GEN)x[3]);
1911: p4=gmul(gneg_i((GEN)p1[2]),p3);
1912:
1913: if (gcmp0((GEN)p1[3]))
1914: {
1915: tetpil=avma;
1916: z[2]=lpile(l,tetpil,gadd(p4,p2));
1917: l=avma; p2=gmul((GEN)x[2],(GEN)x[3]); tetpil=avma;
1918: z[3]=lpile(l,tetpil,gmul2n(p2,1));
1919: copyifstack(p1,z[1]); return z;
1920: }
1921:
1922: p1=gmul((GEN)x[2],(GEN)x[3]);
1923: p1=gmul2n(p1,1); tetpil=avma;
1924: z[2]=ladd(p2,p4); z[3]=ladd(p1,p3);
1925: gerepilemanyvec(l,tetpil,z+2,2);
1926: copyifstack(x[1],z[1]); return z;
1927:
1928: case t_POLMOD:
1929: z=cgetg(lg(x),tx); copyifstack(x[1],z[1]);
1930: l=avma; p1=gsqr((GEN)x[2]); tetpil=avma;
1931: z[2]=lpile(l,tetpil, gres(p1,(GEN)z[1]));
1932: return z;
1933: }
1934:
1935: switch(tx)
1936: {
1937: case t_POL:
1938: {
1939: GEN a = x, p = NULL, pol = NULL;
1940: long vx = varn(x);
1941: av = avma;
1942: if (ff_poltype(&x,&p,&pol))
1943: {
1944: z = quicksqr(x+2, lgef(x)-2);
1945: if (p) z = Fp_pol(z,p);
1946: if (pol) z = from_Kronecker(z,pol);
1947: z = gerepileupto(av, z);
1948: }
1949: else
1950: {
1951: avma = av;
1952: z = quicksqr(a+2, lgef(a)-2);
1953: }
1954: setvarn(z, vx); return z;
1955: }
1956:
1957: case t_SER:
1958: if (gcmp0(x)) return zeroser(varn(x), 2*valp(x));
1959: lx = lg(x); z = cgetg(lx,tx);
1960: z[1] = evalsigne(1) | evalvalp(2*valp(x)) | evalvarn(varn(x));
1961: x += 2; z += 2; lx -= 3;
1962: p2 = (GEN)gpmalloc((lx+1)*sizeof(long));
1963: for (i=0; i<=lx; i++)
1964: {
1965: p2[i] = !isexactzero((GEN)x[i]);
1966: p1=gzero; av=avma; l=(i+1)>>1;
1967: for (j=0; j<l; j++)
1968: if (p2[j] && p2[i-j])
1969: p1 = gadd(p1, gmul((GEN)x[j],(GEN)x[i-j]));
1970: p1 = gshift(p1,1);
1971: if ((i&1) == 0 && p2[i>>1])
1972: p1 = gadd(p1, gsqr((GEN)x[i>>1]));
1973: z[i] = lpileupto(av,p1);
1974: }
1975: z -= 2; free(p2); return normalize(z);
1976:
1977: case t_RFRAC: case t_RFRACN:
1978: z=cgetg(3,tx);
1979: z[1]=lsqr((GEN)x[1]);
1980: z[2]=lsqr((GEN)x[2]); return z;
1981:
1982: case t_QFR: return sqcompreal(x);
1983: case t_QFI: return sqcompimag(x);
1984:
1985: case t_MAT:
1986: lx=lg(x);
1987: if (lx==1) return cgetg(1,tx);
1988: if (lx != lg(x[1])) err(gmuleri,tx,tx);
1989: z=cgetg(lx,tx);
1990: for (j=1; j<lx; j++)
1991: {
1992: z[j]=lgetg(lx,t_COL);
1993: for (i=1; i<lx; i++)
1994: {
1995: p1=gzero; l=avma;
1996: for (k=1; k<lx; k++)
1997: {
1998: p2=gmul(gcoeff(x,i,k),gcoeff(x,k,j));
1999: tetpil=avma; p1=gadd(p1,p2);
2000: }
2001: coeff(z,i,j)=lpile(l,tetpil,p1);
2002: }
2003: }
2004: return z;
2005:
2006: case t_VEC: case t_COL:
2007: err(gmuleri,tx,tx);
2008: }
2009: err(typeer,"gsqr");
2010: return NULL; /* not reached */
2011: }
2012:
2013: /*******************************************************************/
2014: /* */
2015: /* LISTS */
2016: /* */
2017: /*******************************************************************/
2018:
2019: GEN
2020: listcreate(long n)
2021: {
2022: GEN list;
2023:
2024: if (n<0) err(talker,"negative length in listcreate");
2025: list=cgetg(n+2,t_LIST); list[1]=evallgef(2);
2026: return list;
2027: }
2028:
2029: static void
2030: listaffect(GEN list, long index, GEN object)
2031: {
2032: if (index<lgef(list) && isclone(list[index]))
2033: gunclone((GEN)list[index]);
2034: list[index]=lclone(object);
2035: }
2036:
2037: void
2038: listkill(GEN list)
2039: {
2040: long lx,i;
2041:
2042: if (typ(list)!=t_LIST) err(typeer,"listkill");
2043: lx=lgef(list);
2044: for (i=2; i<lx; i++)
2045: if (isclone(list[i])) gunclone((GEN)list[i]);
2046: list[1]=evallgef(2); return;
2047: }
2048:
2049: GEN
2050: listput(GEN list, GEN object, long index)
2051: {
2052: long lx=lgef(list);
2053:
2054: if (typ(list)!=t_LIST) err(typeer,"listput");
2055: if (index < 0) err(talker,"negative index (%ld) in listput", index);
2056: if (!index || index >= lx-1)
2057: {
2058: index = lx-1; lx++;
2059: if (lx > lg(list))
2060: err(talker,"no more room in this list (size %ld)", lg(list)-2);
2061: }
2062: listaffect(list,index+1,object);
2063: list[1]=evallgef(lx);
2064: return (GEN)list[index+1];
2065: }
2066:
2067: GEN
2068: listinsert(GEN list, GEN object, long index)
2069: {
2070: long lx = lgef(list), i;
2071:
2072: if (typ(list)!=t_LIST) err(typeer,"listinsert");
2073: if (index <= 0 || index > lx-1) err(talker,"bad index in listinsert");
2074: lx++; if (lx > lg(list)) err(talker,"no more room in this list");
2075:
2076: for (i=lx-2; i > index; i--) list[i+1]=list[i];
2077: list[index+1] = lclone(object);
2078: list[1] = evallgef(lx); return (GEN)list[index+1];
2079: }
2080:
2081: GEN
2082: gtolist(GEN x)
2083: {
2084: long tx,lx,i;
2085: GEN list;
2086:
2087: if (!x) { list = cgetg(2, t_LIST); list[1] = evallgef(2); return list; }
2088:
2089: tx = typ(x);
2090: lx = (tx==t_LIST)? lgef(x): lg(x);
2091: switch(tx)
2092: {
2093: case t_VEC: case t_COL:
2094: lx++; x--; /* fall through */
2095: case t_LIST:
2096: list = cgetg(lx,t_LIST);
2097: for (i=2; i<lx; i++) list[i] = lclone((GEN)x[i]);
2098: break;
2099: default: err(typeer,"gtolist");
2100: }
2101: list[1] = evallgef(lx); return list;
2102: }
2103:
2104: GEN
2105: listconcat(GEN list1, GEN list2)
2106: {
2107: long i,l1,lx;
2108: GEN list;
2109:
2110: if (typ(list1)!=t_LIST || typ(list2)!=t_LIST)
2111: err(typeer,"listconcat");
2112: l1=lgef(list1)-2; lx=l1+lgef(list2);
2113: list = listcreate(lx-2);
2114: for (i=2; i<=l1+1; i++) listaffect(list,i,(GEN)list1[i]);
2115: for ( ; i<lx; i++) listaffect(list,i,(GEN)list2[i-l1]);
2116: list[1]=evallgef(lx); return list;
2117: }
2118:
2119: GEN
2120: listsort(GEN list, long flag)
2121: {
2122: long i,av=avma, c=list[1], lx = lgef(list)-1;
2123: GEN perm,vec,l;
2124:
2125: if (typ(list) != t_LIST) err(typeer,"listsort");
2126: vec = list+1;
2127: vec[0] = evaltyp(t_VEC) | evallg(lx);
2128: perm = sindexsort(vec);
2129: list[1] = c; l = cgetg(lx,t_VEC);
2130: for (i=1; i<lx; i++) l[i] = vec[perm[i]];
2131: if (flag)
2132: {
2133: c=1; vec[1] = l[1];
2134: for (i=2; i<lx; i++)
2135: if (!gegal((GEN)l[i], (GEN)vec[c]))
2136: vec[++c] = l[i];
2137: else
2138: if (isclone(l[i])) gunclone((GEN)l[i]);
2139: setlgef(list,c+2);
2140: }
2141: else
2142: for (i=1; i<lx; i++) vec[i] = l[i];
2143:
2144: avma=av; return list;
2145: }