Annotation of OpenXM/src/kan96xx/Kan/poly2.c, Revision 1.4
1.4 ! takayama 1: /* $OpenXM: OpenXM/src/kan96xx/Kan/poly2.c,v 1.3 2001/05/04 01:06:25 takayama Exp $ */
1.1 maekawa 2: #include <stdio.h>
3: #include "datatype.h"
4: #include "stackm.h"
5: #include "extern.h"
6: #include "extern2.h"
7:
8: static POLY mapZmonom(POLY f,struct ring *ringp);
9:
10: POLY ppAdd(f,g)
1.3 takayama 11: POLY f; POLY g; /* The result is read only. */
1.1 maekawa 12: {
13: POLY node;
14: struct listPoly nod;
15: POLY h;
16: struct coeff *c;
17: int gt;
18:
19: node = &nod;
20: node->next = POLYNULL;
21: h = node;
22: if (f == POLYNULL) return(g);
23: if (g == POLYNULL) return(f);
24: checkRing(f,g);
25:
26: while (f != POLYNULL && g != POLYNULL) {
27: /*printf("%s + %s\n",POLYToString(f,'*',1),POLYToString(g,'*',1));*/
28: checkRing2(f,g); /* for debug */
29: gt = (*mmLarger)(f,g);
30: switch (gt) {
31: case 1: /* f > g */
32: h -> next = newCell(f->coeffp,f->m);
33: h = h->next;
34: f = f->next;
35: if (f == POLYNULL) {
1.3 takayama 36: h->next = g;
37: return(node->next);
1.1 maekawa 38: }
39: break;
40: case 0: /* f < g */
41: h->next = newCell(g->coeffp,g->m);
42: h = h->next;
43: g = g->next;
44: if (g == POLYNULL) {
1.3 takayama 45: h->next = f;
46: return(node->next);
1.1 maekawa 47: }
48: break;
49: case 2:/* f == g */
50: c = coeffCopy(f->coeffp);
51: Cadd(c,c,g->coeffp);
52: if (!isZero(c)) {
1.3 takayama 53: h->next = newCell(c,f->m);
54: h = h->next;
55: f = f->next;
56: g = g->next;
57: if (f == POLYNULL) {
58: h->next = g;
59: return(node->next);
60: }
61: if (g == POLYNULL) {
62: h->next = f;
63: return(node->next);
64: }
1.1 maekawa 65: }else{
1.3 takayama 66: f = f->next;
67: g = g->next;
68: if (f == POLYNULL) {
69: h->next = g;
70: return(node->next);
71: }
72: if (g == POLYNULL) {
73: h->next = f;
74: return(node->next);
75: }
1.1 maekawa 76: }
77: break;
78: default:
79: errorPoly("ppAdd(). Internal error. Invalid value of mmLarger.");
80: break;
81: }
82: }
83: return(node->next);
84: }
85:
86: POLY ppSub(f,g)
1.3 takayama 87: POLY f; POLY g; /* The result is read only. */
1.1 maekawa 88: {
89: POLY h;
90: struct coeff *c;
91:
92: if (g == POLYNULL) return(f);
93: if (f == POLYNULL) return(cpMult(intToCoeff(-1,g->m->ringp),g));
94: checkRing(f,g);
95:
96: h = cpMult(intToCoeff(-1,g->m->ringp),g);
97: return(ppAdd(f,h));
98: }
99:
100:
101: POLY cpMult(c,f)
1.3 takayama 102: struct coeff *c;
103: POLY f;
1.1 maekawa 104: {
105: POLY node;
106: struct listPoly nod;
107: POLY h;
108: struct coeff *newc;
109: int p;
110: node = &nod;
111: if (f == POLYNULL || isZero(c)) return(POLYNULL);
112: p = f->m->ringp->p;
113: node ->next = POLYNULL;
114: h = node;
115: while (f != POLYNULL) {
116: newc = coeffCopy(c);
117: Cmult(newc,newc,f->coeffp);
118: if ((p==0) || !isZero(newc)) {
119: h->next = newCell(newc,f->m);
120: h = h->next;
121: }
122: f = f->next;
123: }
124: return(node->next);
125: }
126:
127: MONOMIAL monomialAdd_poly(m,m2)
1.3 takayama 128: MONOMIAL m,m2;
1.1 maekawa 129: {
130: extern int Msize;
131: MONOMIAL f;
132: int i;
133: int n;
134: n = m->ringp->n;
135: f = (MONOMIAL) sGC_malloc(sizeof(struct smallMonomial)+n*Msize);
136:
137: if (f == (MONOMIAL) NULL) errorPoly("No more memory.");
138: f->ringp = m->ringp;
139: for (i=0; i<n; i++) {
140: (f->e)[i].x = (m->e)[i].x + (m2->e)[i].x;
141: (f->e)[i].D = (m->e)[i].D + (m2->e)[i].D;
142: }
143: return(f);
144: }
145:
146: /* Note that mpMult_poly is called from mmLarger_tower! */
147: POLY mpMult_poly(f,g)
1.3 takayama 148: POLY f;
149: POLY g;
1.1 maekawa 150: {
151: POLY node;
152: struct listPoly nod;
153: struct coeff *c;
154: int n,i;
155: POLY h;
156: MONOMIAL m;
157: int p;
158: node = &nod;
159: if (f == POLYNULL || g == POLYNULL) return(POLYNULL);
160: node->next = POLYNULL;
161: h = node;
162: checkRing(f,g);
163: n = f->m->ringp->n; p = f->m->ringp->p;
164: while(g != POLYNULL) {
165: checkRing2(f,g);
166: c = coeffCopy(f->coeffp);
167: Cmult(c,c,g->coeffp);
168: if ((p==0) || !isZero(c)) {
169: m = (*monomialAdd)(f->m,g->m);
170: h->next = newCell(c,m);
171: h = h->next;
172: }
173: g = g->next;
174: }
175: return(node->next);
176: }
177:
178: POLY ppMult_old(f,g)
1.3 takayama 179: POLY f,g;
1.1 maekawa 180: {
181: POLY r;
182: POLY tmp;
183: r = POLYNULL;
184: while( f != POLYNULL) {
185: tmp = (*mpMult)(f,g);
186: r = ppAddv(r,tmp); /* r and tmp will be broken */
187: f = f->next;
188: }
189: return(r);
190: }
191:
192: POLY ppAddv(f,g)
1.3 takayama 193: POLY f; POLY g; /* It breaks f and g. Use it just after calling mpMult() */
1.1 maekawa 194: {
195: POLY node;
196: struct listPoly nod;
197: POLY h;
198: struct coeff *c;
199: int gt;
200:
201: node = &nod;
202: /*printf("ppAddv:1%s + %s\n",POLYToString(f,'*',1),POLYToString(g,'*',1));*/
203: node->next = POLYNULL;
204: h = node;
205: if (f == POLYNULL) return(g);
206: if (g == POLYNULL) return(f);
207: checkRing(f,g);
208:
209: while (f != POLYNULL && g != POLYNULL) {
210: checkRing2(f,g); /* for debug */
211: /*printf("%s + %s\n",POLYToString(f,'*',1),POLYToString(g,'*',1));*/
212: gt = (*mmLarger)(f,g);
213: switch (gt) {
214: case 1: /* f > g */
215: h->next = f;
216: h = h->next; f = f->next;;
217: if (f == POLYNULL) {
1.3 takayama 218: h->next = g;
219: return(node->next);
1.1 maekawa 220: }
221: break;
222: case 0: /* f < g */
223: h->next = g;
224: h = h->next; g = g->next;
225: if (g == POLYNULL) {
1.3 takayama 226: h->next = f;
227: return(node->next);
1.1 maekawa 228: }
229: break;
230: case 2:/* f == g */
231: c = f->coeffp;
232: Cadd(c,c,g->coeffp);
233: if (!isZero(c)) {
1.3 takayama 234: h->next = f;
235: h = h->next; f = f->next;;
236: g = g->next;
237: if (f == POLYNULL) {
238: h->next = g;
239: return(node->next);
240: }
241: if (g == POLYNULL) {
242: h->next = f;
243: return(node->next);
244: }
1.1 maekawa 245: }else{
1.3 takayama 246: f = f->next;
247: g = g->next;
248: if (f == POLYNULL) {
249: h->next = g;
250: return(node->next);
251: }
252: if (g == POLYNULL) {
253: h->next = f;
254: return(node->next);
255: }
1.1 maekawa 256: }
257: break;
258: default:
259: errorPoly("ppAddv(). Internal error. Invalid value of mmLarger.");
260: break;
261: }
262: }
263: return(node->next);
264: }
265:
266: POLY pPower(f,k)
1.3 takayama 267: POLY f;
268: int k;
1.1 maekawa 269: {
270: POLY r;
271: int i,n;
272: if (f == POLYNULL) return(POLYNULL); /* Is it ok? 0^0 = 0.*/
273: if (k == 0) return(cxx(1,0,0,f->m->ringp));
274: if (f->next == POLYNULL && k<0) {
275: /* when f is monomial. */
276: r = newCell(coeffCopy(f->coeffp),monomialCopy(f->m));
277: n = r->m->ringp->n;
278: for (i=0; i<n; i++) {
279: r->m->e[i].x *= k;
280: r->m->e[i].D *= k;
281: }
282: if (!isOne(r->coeffp)) {
283: warningPoly("pPower(poly,negative integer) not implemented yet. Returns 1.");
284: r = cxx(1,0,0,f->m->ringp);
285: }
286: return(r);
287: }
288: r = cxx(1,0,0,f->m->ringp);
289: if (k < 0) {
290: warningPoly("pPower(poly,negative integer) not implemented yet. Returns 1.");
291: }
292: for (i=0; i<k; i++) {
293: r = ppMult(f,r);
294: }
295: return(r);
296: }
297:
298: POLY pPower_poly(f,k)
1.3 takayama 299: POLY f;
300: int k;
1.1 maekawa 301: {
302: POLY r;
303: int i,n;
304: if (f == POLYNULL) return(POLYNULL); /* Is it ok? 0^0 = 0.*/
305: if (k == 0) return(cxx(1,0,0,f->m->ringp));
306: if (f->next == POLYNULL && k<0) {
307: /* when f is monomial. */
308: r = newCell(coeffCopy(f->coeffp),monomialCopy(f->m));
309: n = r->m->ringp->n;
310: for (i=0; i<n; i++) {
311: r->m->e[i].x *= k;
312: r->m->e[i].D *= k;
313: }
314: if (!isOne(r->coeffp)) {
315: warningPoly("pPower_poly(poly,negative integer) not implemented yet. Returns 1.");
316: r = cxx(1,0,0,f->m->ringp);
317: }
318: return(r);
319: }
320: r = cxx(1,0,0,f->m->ringp);
321: if (k < 0) {
322: warningPoly("pPower_poly(poly,negative integer) not implemented yet. Returns 1.");
323: }
324: for (i=0; i<k; i++) {
325: r = ppMult_poly(f,r);
326: }
327: return(r);
328: }
329:
330: POLY modulop_trash(f,ringp)
1.3 takayama 331: POLY f;
332: struct ring *ringp;
1.1 maekawa 333: {
334: int p;
335: POLY h;
336: MP_INT *c;
337: int cc;
338: POLY node;
339: struct ring *nextRing;
340: POLY fc;
341:
342: if (f == POLYNULL) return(f);
343: p = ringp->p;
344: if (f->m->ringp->p != 0) {
345: warningPoly("modulop(f,ringp) must be called with f in the characteristic 0 ring. Returns 0.");
346: return(POLYNULL);
347: }
348: if (f->m->ringp->n != ringp->n) {
349: warningPoly("modulop(f,ringp): f->m->ringp->n must be equal to ringp->n. Returns 0.");
350: return(POLYNULL);
351: }
352:
353: /* The case of ringp->next != NULL */
354: if (ringp->next != (struct ring *)NULL) {
355: nextRing = ringp->next;
356: node = newCell(newCoeff(),newMonomial(ringp));
357: node->next = POLYNULL;
358: h = node;
359:
360: while (f != POLYNULL) {
361: fc = bxx(f->coeffp->val.bigp,0,0,nextRing);
362: h->next = newCell(newCoeff(),monomialCopy(f->m));
363: h = h->next;
364: h->m->ringp = ringp;
365: h->coeffp->tag = POLY_COEFF;
366: h->coeffp->p = p;
367: h->coeffp->val.f = fc;
368: f = f->next;
369: }
370: return(node->next);
371: }
372:
373:
374: /* In case of ringp->next == NULL */
375: if (p != 0) {
376: c = newMP_INT();
377: node = newCell(newCoeff(),newMonomial(ringp));
378: node->next = POLYNULL;
379: h = node;
380:
381: while (f != POLYNULL) {
382: mpz_mod_ui(c,f->coeffp->val.bigp,(unsigned long int)p);
383: cc = (int) mpz_get_si(c);
384: if (cc != 0) {
1.3 takayama 385: h->next = newCell(newCoeff(),monomialCopy(f->m));
386: h = h->next;
387: h->m->ringp = ringp;
388: h->coeffp->tag = INTEGER;
389: h->coeffp->p = p;
390: h->coeffp->val.i = cc;
1.1 maekawa 391: }
392: f = f->next;
393: }
394: return(node->next);
395: }else{
396: h = f = pcmCopy(f);
397: while (f != POLYNULL) {
398: f->m->ringp = ringp;
399: f = f->next;
400: }
401: return(h);
402: }
403:
404: }
405:
406: POLY modulop(f,ringp)
1.3 takayama 407: POLY f;
408: struct ring *ringp;
409: /* Z[x] ---> R[x] where R=Z, Z/Zp, ringp->next. */
1.1 maekawa 410: {
411: int p;
412: POLY h;
413: MP_INT *c;
414: int cc;
415: POLY node;
416: POLY fc;
417:
418: if (f == POLYNULL) return(f);
419: p = ringp->p;
420: if (f->m->ringp->p != 0 || f->m->ringp->next != (struct ring *)NULL) {
421: warningPoly("modulop(f,ringp) must be called with f in the characteristic 0 ring Z[x]. Returns 0.");
422: return(POLYNULL);
423: }
424: if (f->m->ringp->n != ringp->n) {
425: warningPoly("modulop(f,ringp): f->m->ringp->n must be equal to ringp->n. Returns 0.");
426: return(POLYNULL);
427: }
428:
429: /* [1] The case of R = ringp->next */
430: if (ringp->next != (struct ring *)NULL) {
431: h = ZERO;
432: while (f != POLYNULL) {
433: h = ppAdd(h,mapZmonom(f,ringp));
434: f = f->next;
435: }
436: return(h);
437: }
438:
439: /* [2] The case of R = Z/Zp */
440: if (p != 0) {
441: c = newMP_INT();
442: node = newCell(newCoeff(),newMonomial(ringp));
443: node->next = POLYNULL;
444: h = node;
445:
446: while (f != POLYNULL) {
447: mpz_mod_ui(c,f->coeffp->val.bigp,(unsigned long int)p);
448: cc = (int) mpz_get_si(c);
449: if (cc != 0) {
1.3 takayama 450: h->next = newCell(newCoeff(),monomialCopy(f->m));
451: h = h->next;
452: h->m->ringp = ringp;
453: h->coeffp->tag = INTEGER;
454: h->coeffp->p = p;
455: h->coeffp->val.i = cc;
1.1 maekawa 456: }
457: f = f->next;
458: }
459: return(node->next);
460: }
461:
462: /* [3] The case of R = Z */
463: h = f = pcmCopy(f);
464: while (f != POLYNULL) {
465: f->m->ringp = ringp;
466: f = f->next;
467: }
468: return(h);
469:
470:
471: }
472:
473: POLY modulopZ(f,pcoeff)
1.3 takayama 474: POLY f;
475: struct coeff *pcoeff;
476: /* Z[x] ---> Z[x] , f ---> f mod pcoeff*/
1.1 maekawa 477: {
478: int p;
479: POLY h;
480: struct coeff *c;
481: int cc;
482: POLY node;
483: POLY fc;
484: MP_INT *bigp;
485: MP_INT *tmp;
486: struct ring *ringp;
487:
488: if (f == POLYNULL) return(f);
489: ringp = f->m->ringp;
490: if (pcoeff->tag != MP_INTEGER) {
491: warningPoly("modulopZ(): pcoeff must be a universalNumber.");
492: warningPoly("Returns 0.");
493: return(POLYNULL);
494: }
495: bigp = pcoeff->val.bigp;
496: if (f->m->ringp->p != 0 || f->m->ringp->next != (struct ring *)NULL) {
497: warningPoly("modulopZ(f,p) must be called with f in the characteristic 0 ring Z[x]. Returns 0.");
498: return(POLYNULL);
499: }
500: if (f->m->ringp->n != ringp->n) {
501: warningPoly("modulop(f,ringp): f->m->ringp->n must be equal to ringp->n. Returns 0.");
502: return(POLYNULL);
503: }
504:
505: /* [1] The case of R = ringp->next */
506: if (ringp->next != (struct ring *)NULL) {
507: warningPoly("modulopZ workds only for flat polynomials. Returns 0.");
508: return(POLYNULL);
509: }
510:
511: /* [2] The case of R = Z */
512: node = newCell(newCoeff(),newMonomial(ringp));
513: node->next = POLYNULL;
514: h = node;
515:
516: c = newCoeff();
517: tmp = newMP_INT();
518: while (f != POLYNULL) {
519: mpz_mod(tmp,f->coeffp->val.bigp,bigp);
520: if (mpz_sgn(tmp) != 0) {
521: c->tag = MP_INTEGER;
522: c->p = 0;
523: c->val.bigp = tmp;
524: h->next = newCell(c,monomialCopy(f->m));
525: h = h->next;
526: h->m->ringp = ringp;
527:
528: c = newCoeff();
529: tmp = newMP_INT();
530: }
531: f = f->next;
532: }
533: return(node->next);
534:
535: }
536:
537: struct pairOfPOLY quotientByNumber(f,pcoeff)
1.3 takayama 538: POLY f;
539: struct coeff *pcoeff;
540: /* Z[x] ---> Z[x],Z[x] , f = first*pcoeff + second */
1.1 maekawa 541: {
542: int p;
543: POLY h;
544: POLY h2;
545: struct coeff *c;
546: int cc;
547: POLY node;
548: struct coeff *c2;
549: POLY node2;
550: POLY fc;
551: MP_INT *bigp;
552: MP_INT *tmp;
553: MP_INT *tmp2;
554: struct ring *ringp;
555: struct pairOfPOLY r;
556:
557: if (f == POLYNULL) {
558: r.first = f; r.second = f;
559: return(r);
560: }
561: ringp = f->m->ringp;
562: if (pcoeff->tag != MP_INTEGER) {
563: warningPoly("quotientByNumber(): pcoeff must be a universalNumber.");
564: warningPoly("Returns (0,0).");
565: r.first = f; r.second = f;
566: return(r);
567: }
568: bigp = pcoeff->val.bigp;
569: if (f->m->ringp->p != 0 || f->m->ringp->next != (struct ring *)NULL) {
570: warningPoly("quotientByNumber(f,p) must be called with f in the characteristic 0 ring Z[x]. Returns (0,0).");
571: r.first = f; r.second = f;
572: return(r);
573: }
574: if (f->m->ringp->n != ringp->n) {
575: warningPoly("quotientByNumber(f,p): f->m->ringp->n must be equal to ringp->n. Returns 0.");
576: r.first = f; r.second = f;
577: return(r);
578: }
579:
580: /* [1] The case of R = ringp->next */
581: if (ringp->next != (struct ring *)NULL) {
582: warningPoly("quotientByNumber() workds only for flat polynomials. Returns 0.");
583: r.first = f; r.second = f;
584: return(r);
585: }
586:
587: /* [2] The case of R = Z */
588: node = newCell(newCoeff(),newMonomial(ringp));
589: node->next = POLYNULL;
590: h = node;
591: node2 = newCell(newCoeff(),newMonomial(ringp));
592: node2->next = POLYNULL;
593: h2 = node2;
594:
595: c = newCoeff();
596: tmp = newMP_INT();
597: c2 = newCoeff();
598: tmp2 = newMP_INT();
599: while (f != POLYNULL) {
600: mpz_mod(tmp,f->coeffp->val.bigp,bigp);
601: if (mpz_sgn(tmp) != 0) {
602: c->tag = MP_INTEGER;
603: c->p = 0;
604: c->val.bigp = tmp;
605: h->next = newCell(c,monomialCopy(f->m));
606: h = h->next;
607: h->m->ringp = ringp;
608:
609: c = newCoeff();
610: tmp = newMP_INT();
611: }
612: mpz_tdiv_q(tmp2,f->coeffp->val.bigp,bigp);
613: if (mpz_sgn(tmp2) != 0) {
614: c2->tag = MP_INTEGER;
615: c2->p = 0;
616: c2->val.bigp = tmp2;
617: h2->next = newCell(c2,monomialCopy(f->m));
618: h2 = h2->next;
619: h2->m->ringp = ringp;
620:
621: c2 = newCoeff();
622: tmp2 = newMP_INT();
623: }
624: f = f->next;
625: }
626: r.first = node2->next;
627: r.second = node->next;
628: return(r);
629:
630: }
631:
632:
633: POLY modulo0(f,ringp)
1.3 takayama 634: POLY f;
635: struct ring *ringp;
1.1 maekawa 636: {
637: int p;
638: POLY h;
639: struct coeff *c;
640: POLY node;
641: if (f == POLYNULL) return(f);
642: p = ringp->p;
643: if (p != 0) {
644: warningPoly("modulo0(f,ringp) must be called with the characteristic 0 ring*ringp. Returns 0.");
645: return(POLYNULL);
646: }
647: switch (f->coeffp->tag) {
648: case MP_INTEGER:
649: if (f->m->ringp->p == 0) {
650: node = pcmCopy(f);
651: f = node;
652: while (f != POLYNULL) {
1.3 takayama 653: f->m->ringp = ringp; /* Touch the monomial "ringp" field. */
654: f = f->next;
1.1 maekawa 655: }
656: return(node);
657: }
658: break;
659: case POLY_COEFF:
660: node = pcmCopy(f);
661: f = node;
662: while (f != POLYNULL) {
663: f->m->ringp = ringp; /* Touch the monomial "ringp" field. */
664: f = f->next;
665: }
666: return(node);
667: break;
668: case INTEGER:
669: node = newCell(newCoeff(),newMonomial(ringp));
670: node->next = POLYNULL;
671: h = node;
672:
673: while (f != POLYNULL) {
674: c = newCoeff();
675: c->tag = MP_INTEGER;
676: c->p = 0;
677: c->val.bigp = newMP_INT();
678: mpz_set_si(c->val.bigp,f->coeffp->val.i);
679: h->next = newCell(c,monomialCopy(f->m));
680: h = h->next;
681: h->m->ringp = ringp;
682: f = f->next;
683: }
684: return(node->next);
685: break;
686: default:
687: warningPoly("modulo0(): coefficients have to be MP_INTEGER or INTEGER. Returns 0");
688: return(POLYNULL);
689: break;
690: }
691: }
692:
693:
694: struct object test(ob) /* test3 */
1.3 takayama 695: struct object ob;
1.1 maekawa 696: {
697: struct object rob;
698: int k;
699: static POLY f0;
700: static POLY f1;
701:
702: POLY addNode;
703: POLY f;
704: int i;
705: static int s=0;
706: MP_INT *mp;
707: extern struct ring *SmallRingp;
708: extern struct ring *CurrentRingp;
709: addNode = pMalloc(SmallRingp);
710: k = ob.lc.ival;
711: switch(s) {
712: case 0:
713: f0 = addNode;
714: for (i=k; i>=0; i--) {
715: f0->next = bxx(BiiPower(-k,i),0,i,CurrentRingp);
716: if (f0->next != POLYNULL) {
1.3 takayama 717: f0 = f0->next;
1.1 maekawa 718: }
719: }
720: f0 = addNode->next;
721: s++;
722: rob.lc.poly = f0;
723: break;
724: case 1:
725: f1 = addNode;
726: for (i=k; i>=0; i--) {
727: f1->next = bxx(BiiPower(k,i),0,i,CurrentRingp);
728: if (f1->next != POLYNULL) {
1.3 takayama 729: f1 = f1->next;
1.1 maekawa 730: }
731: }
732: f1 = addNode->next;
733: s = 0;
734: rob.lc.poly = f1;
735: break;
736: default:
737: rob.lc.poly = POLYNULL;
738: s = 0;
739: break;
740: }
741:
742:
743: rob.tag = Spoly;
744: return(rob);
745: }
746:
747:
748: int pLength(f)
1.3 takayama 749: POLY f;
1.1 maekawa 750: {
751: int c=0;
752: if (f ISZERO) return(0);
753: while (f != POLYNULL) {
754: c++;
755: f = f->next;
756: }
757: return(c);
758: }
759:
760:
761: POLY ppAddv2(f,g,top,nexttop)
1.3 takayama 762: POLY f; POLY g; /* It breaks f and g. Use it just after calling mpMult() */
763: POLY top;
764: POLY *nexttop;
765: /* top is the starting address in the list f.
766: if top == POLYNULL, start from f.
767:
768: *nexttop == 0
769: == g
770: == h or 0
1.1 maekawa 771:
1.3 takayama 772: It must be called as r = ppAddv2(r,g,...);
1.1 maekawa 773: */
774: {
775: POLY node;
776: struct listPoly nod;
777: POLY h;
778: struct coeff *c;
779: int gt;
780: POLY g0;
781: POLY f0; /* for debug */
782:
783: node = &nod;
784: /* printf("ppAddv:1%s + %s\n",POLYToString(f,'*',1),POLYToString(g,'*',1)); */
785: node->next = POLYNULL;
786: h = node;
787: *nexttop = POLYNULL;
788: if (f == POLYNULL) return(g);
789: if (g == POLYNULL) return(f);
790: checkRing(f,g);
791:
792: f0 = f;
793: if (top != POLYNULL) {
794: while (f != top) {
795: if (f == POLYNULL) {
1.3 takayama 796: fprintf(stderr,"\nppAddv2(): Internal error.\n");fflush(stderr);
797: fprintf(stderr,"f = %s\n",POLYToString(f0,'*',0));
798: fprintf(stderr,"g = %s\n",POLYToString(g0,'*',0));
799: fprintf(stderr,"top=%s\n",POLYToString(top,'*',0));
800: errorPoly("ppAddv2(). Internal error=1.");
1.1 maekawa 801: }
802: h->next = f;
803: h = h->next;
804: f = f->next;
805: }
806: }
807: g0 = g;
808: *nexttop = g0;
809:
810: while (f != POLYNULL && g != POLYNULL) {
811: checkRing2(f,g); /* for debug */
812: /* printf("%s + %s\n",POLYToString(f,'*',1),POLYToString(g,'*',1)); */
813: gt = (*mmLarger)(f,g);
814: switch (gt) {
815: case 1: /* f > g */
816: h->next = f;
817: h = h->next; f = f->next;;
818: if (f == POLYNULL) {
1.3 takayama 819: h->next = g;
820: return(node->next);
1.1 maekawa 821: }
822: break;
823: case 0: /* f < g */
824: h->next = g;
825: h = h->next; g = g->next;
826: if (g == POLYNULL) {
1.3 takayama 827: h->next = f;
828: return(node->next);
1.1 maekawa 829: }
830: break;
831: case 2:/* f == g */
832: c = g->coeffp;
833: Cadd(c,f->coeffp,c);
834: if (!isZero(c)) {
1.3 takayama 835: h->next = g;
836: h = h->next;
837: f = f->next;;
838: g = g->next;
839: if (f == POLYNULL) {
840: h->next = g;
841: return(node->next);
842: }
843: if (g == POLYNULL) {
844: h->next = f;
845: return(node->next);
846: }
1.1 maekawa 847: }else{
1.3 takayama 848: if (g == g0) {
849: if (h != node) {
850: *nexttop = h;
851: }else{
852: *nexttop = POLYNULL;
853: }
854: }
1.1 maekawa 855:
1.3 takayama 856: f = f->next;
857: g = g->next;
1.1 maekawa 858:
1.3 takayama 859: if (f == POLYNULL) {
860: h->next = g;
861: return(node->next);
862: }
863: if (g == POLYNULL) {
864: h->next = f;
865: return(node->next);
866: }
1.1 maekawa 867: }
868: break;
869: default:
870: errorPoly("ppAddv(). Internal error. Invalid value of mmLarger.");
871: break;
872: }
873: }
874: return(node->next);
875: }
876:
877: POLY ppMult(f,g)
1.3 takayama 878: POLY f,g;
1.1 maekawa 879: {
880: POLY r;
881: POLY tmp;
882: POLY top;
883: POLY nexttop;
884: r = POLYNULL; top = POLYNULL;
885: while( f != POLYNULL) {
886: /* tmp = (*mpMult)(f,g); (*mpMult) is no more used. */
887: tmp = (*(f->m->ringp->multiplication))(f,g);
888: /*printf("mpMult(%s,%s) ->%s\n",POLYToString(f,'*',1),POLYToString(g,'*',1),POLYToString(tmp,'*',1)); */
889: r = ppAddv2(r,tmp,top,&nexttop); /* r and tmp will be broken */
890: top = nexttop;
891: f = f->next;
892: }
893: return(r);
894: }
895:
896: POLY ppMult_poly(f,g)
1.3 takayama 897: POLY f,g;
1.1 maekawa 898: {
899: POLY r;
900: POLY tmp;
901: POLY top;
902: POLY nexttop;
903: r = POLYNULL; top = POLYNULL;
904: while( f != POLYNULL) {
905: tmp = mpMult_poly(f,g);
906: r = ppAddv2(r,tmp,top,&nexttop); /* r and tmp will be broken */
907: top = nexttop;
908: f = f->next;
909: }
910: return(r);
911: }
912:
913: POLY mapZmonom(f,ringp)
1.3 takayama 914: POLY f; /* assumes monomial. f \in Z[x] */
915: struct ring *ringp; /* R[x] */
1.1 maekawa 916: {
917: struct ring *nextRing;
918: struct ring nextRing0;
919: POLY ff;
920: POLY node;
921: POLY gg;
922: int l,c,d;
923:
924: nextRing = ringp->next;
925: nextRing0 = *nextRing; nextRing0.p = 0; nextRing0.next = 0;
926: /* nextRing0 = Z[y] where y is the variables of R. */
927:
928: ff = bxx(f->coeffp->val.bigp,0,0,&nextRing0);
929: ff = modulop(ff,nextRing);
930:
931: node = newCell(newCoeff(),monomialCopy(f->m));
932: node->next = POLYNULL;
933: node->m->ringp = ringp;
934:
935: node->coeffp->p = ringp->p;
936:
937: l = ringp->l; c = ringp->c;
938: /* If q-analog q x ---> (q) x. */
939: if (l-c > 0) {
940: d = node->m->e[0].x; /* degree of q in R[x].*/
941: node->m->e[0].x = 0;
942: gg = cxx(1,0,d,nextRing); /* q^d = x[0]^d */
943: ff = ppMult(gg,ff);
944: }
945:
946: node->coeffp->tag = POLY_COEFF;
947: node->coeffp->val.f = ff;
948: return(node);
949: }
1.4 ! takayama 950:
! 951: POLY reduceContentOfPoly(POLY f,struct coeff **contp) {
! 952: struct coeff *cont;
! 953: struct coeff *cOne = NULL;
! 954: extern struct ring *SmallRingp;
! 955: if (cOne == NULL) cOne = intToCoeff(1,SmallRingp);
! 956:
! 957: if (f == POLYNULL) return f;
! 958: if (f->m->ringp->p != 0) return f;
! 959: if (f->coeffp->tag != MP_INTEGER) return f;
! 960: cont = gcdOfCoeff(f);
! 961: *contp = cont;
! 962: if (coeffGreater(cont,cOne)) {
! 963: f = quotientByNumber(f,cont).first;
! 964: }
! 965: return f;
! 966: }
! 967:
! 968: int coeffSizeMin(POLY f) {
! 969: int size;
! 970: int t;
! 971: if (f == POLYNULL) return 0;
! 972: if (f->m->ringp->p != 0) return 0;
! 973: if (f->coeffp->tag != MP_INTEGER) return 0;
! 974: size = mpz_size(f->coeffp->val.bigp);
! 975: while (f != POLYNULL) {
! 976: t = mpz_size(f->coeffp->val.bigp);
! 977: if (t < size) size = t;
! 978: if (size == 1) return size;
! 979: f = f->next;
! 980: }
! 981: }
! 982:
! 983: struct coeff *gcdOfCoeff(POLY f) {
! 984: extern struct ring *SmallRingp;
! 985: struct coeff *t;
! 986: MP_INT *tmp;
! 987: MP_INT *tmp2;
! 988: static MP_INT *cOne = NULL;
! 989: if (cOne == NULL) {
! 990: cOne = newMP_INT();
! 991: mpz_set_si(cOne,(long) 1);
! 992: }
! 993: if (f == POLYNULL) return 0;
! 994: if (f->m->ringp->p != 0) return intToCoeff(0,SmallRingp);
! 995: if (f->coeffp->tag != MP_INTEGER) return intToCoeff(0,SmallRingp);
! 996: tmp = f->coeffp->val.bigp;
! 997: tmp2 = newMP_INT();
! 998: while (f != POLYNULL) {
! 999: mpz_gcd(tmp2,tmp,f->coeffp->val.bigp); /* tmp = tmp2 OK? */
! 1000: tmp = tmp2;
! 1001: if (mpz_cmp(tmp,cOne)==0) return intToCoeff(1,SmallRingp);
! 1002: f = f->next;
! 1003: }
! 1004: return mpintToCoeff(tmp,SmallRingp);
1.1 maekawa 1005:
1.4 ! takayama 1006: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to poly2.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM / src / kan96xx / Kan