Annotation of OpenXM/src/kan96xx/Kan/order.c, Revision 1.17
1.17 ! takayama 1: /* $OpenXM: OpenXM/src/kan96xx/Kan/order.c,v 1.16 2018/09/07 00:15:44 takayama Exp $ */
1.1 maekawa 2: #include <stdio.h>
1.15 ohara 3: #include <stdlib.h>
1.17 ! takayama 4: #include <string.h>
1.1 maekawa 5: #include "datatype.h"
6: #include "stackm.h"
7: #include "extern.h"
8: #include "extern2.h"
9:
10: /* The format of order.
11: Example: graded lexicographic order
12: x_{N-1} x_{N-2} ... x_0 D_{N-1} .... D_{0}
13: 1 1 1 1 1
14: 1 0 0 0 0
15: 0 1 0 0 0
16: ..............................................
17:
18: (ringp->order)[i][j] should be (ringp->order)[i*2*N+j].
19: All order matrix is generated by functions in smacro.sm1
20: */
21:
22: static void warningOrder(char *s);
23: static void errorOrder(char *s);
24:
25: void setOrderByMatrix(order,n,c,l,omsize)
1.4 takayama 26: int order[];
27: int n,c,l,omsize;
1.1 maekawa 28: {
29: int i,j;
30: int *Order;
31: extern struct ring *CurrentRingp;
32:
33: switch_mmLarger("default");
1.4 takayama 34: /* q-case */
1.1 maekawa 35: if ( l-c > 0) {
36: switch_mmLarger("qmatrix");
37: }
38:
39: Order = (int *)sGC_malloc(sizeof(int)*(2*n)*(omsize));
40: if (Order == (int *)NULL) errorOrder("No memory.");
41: CurrentRingp->order = Order;
42: CurrentRingp->orderMatrixSize = omsize;
43: for (i=0; i<omsize; i++) {
44: for (j=0; j<2*n; j++) {
45: Order[i*2*n+j] = order[i*2*n+j];
46: }
47: }
48: }
49:
50: void showRing(level,ringp)
1.4 takayama 51: int level;
52: struct ring *ringp;
1.1 maekawa 53: {
54: int i,j;
55: FILE *fp;
56: char tmp[100];
57: int N,M,L,C,NN,MM,LL,CC;
58: char **TransX,**TransD;
59: int *Order;
60: int P;
61: char *mtype;
62: extern char *F_isSameComponent;
1.5 takayama 63: POLY f;
1.6 takayama 64: POLY fx;
65: POLY fd;
66: POLY rf;
1.1 maekawa 67: fp = stdout;
68:
69: N=ringp->n; M = ringp->m; L = ringp->l; C = ringp->c;
70: NN=ringp->nn; MM = ringp->mm; LL = ringp->ll; CC = ringp->cc;
71: TransX = ringp->x; TransD = ringp->D;
72: Order = ringp->order;
73: P = ringp->p;
74:
75:
76: fprintf(fp,"\n---------- the current ring ---- name: %s------\n",ringp->name);
77: fprintf(fp,"Characteristic is %d. ",P);
78: fprintf(fp,"N0=%d N=%d NN=%d M=%d MM=%d L=%d LL=%d C=%d CC=%d omsize=%d\n",N0,N,NN,M,MM,L,LL,C,CC,ringp->orderMatrixSize);
79: fprintf(fp,"\n");
80:
81: /* print identifier names */
82: if (N-M >0) {
83: fprintf(fp,"Differential variables: ");
84: for (i=M; i<N; i++) fprintf(fp," %4s ",TransX[i]);
85: for (i=M; i<N; i++) fprintf(fp," %4s ",TransD[i]);
86: fprintf(fp,"\n");
87: fprintf(fp,"where ");
88: for (i=M; i<N; i++) {
1.6 takayama 89: fx = cxx(1,i,1,ringp); fd = cdd(1,i,1,ringp);
90: rf = ppSub(ppMult(fd,fx),ppMult(fx,fd));
91: fprintf(fp," %s %s - %s %s = %s, ",TransD[i],TransX[i],
92: TransX[i],TransD[i],POLYToString(rf,'*',0));
1.1 maekawa 93: }
94: fprintf(fp,"\n\n");
95: }
96: if (M-L >0) {
97: fprintf(fp,"Difference variables: ");
98: for (i=L; i<M; i++) fprintf(fp," %4s ",TransX[i]);
99: for (i=L; i<M; i++) fprintf(fp," %4s ",TransD[i]);
100: fprintf(fp,"\n");
101: fprintf(fp,"where ");
102: for (i=L; i<M; i++) {
1.5 takayama 103: fprintf(fp," %s %s - %s %s = ",TransD[i],TransX[i],
104: TransX[i],TransD[i]);
105: f=ppSub(ppMult(cdd(1,i,1,ringp),cxx(1,i,1,ringp)),
106: ppMult(cxx(1,i,1,ringp),cdd(1,i,1,ringp)));
107: fprintf(fp," %s, ",POLYToString(f,'*',0));
1.1 maekawa 108: }
109: fprintf(fp,"\n\n");
110: }
111: if (L-C >0) {
112: fprintf(fp,"q-Difference variables: ");
113: for (i=C; i<L; i++) fprintf(fp," %4s ",TransX[i]);
114: for (i=C; i<L; i++) fprintf(fp," %4s ",TransD[i]);
115: fprintf(fp,"\n");
116: fprintf(fp,"where ");
117: for (i=C; i<L; i++) {
118: fprintf(fp," %s %s = %s %s %s, ",TransD[i],TransX[i],
1.4 takayama 119: TransX[0],
120: TransX[i],TransD[i]);
1.1 maekawa 121: }
122: fprintf(fp,"\n\n");
123: }
124: if (C>0) {
125: fprintf(fp,"Commutative variables: ");
126: for (i=0; i<C; i++) fprintf(fp," %4s ",TransX[i]);
127: for (i=0; i<C; i++) fprintf(fp," %4s ",TransD[i]);
128: fprintf(fp,"\n\n");
129: }
130:
131: if (strcmp(F_isSameComponent,"x") == 0) {
132: fprintf(fp,"Integral or summation or graduation variables are : ");
133: for (i=CC; i<C; i++) fprintf(fp," %4s ",TransX[i]);
134: for (i=LL; i<L; i++) fprintf(fp," %4s ",TransX[i]);
135: for (i=MM; i<M; i++) fprintf(fp," %4s ",TransX[i]);
136: for (i=NN; i<N; i++) fprintf(fp," %4s ",TransX[i]);
137: fprintf(fp,"\n");
138: }else if (strcmp(F_isSameComponent,"xd") == 0) {
139: fprintf(fp,"Graduation variables are : ");
140: for (i=CC; i<C; i++) fprintf(fp," %4s ",TransX[i]);
141: for (i=LL; i<L; i++) fprintf(fp," %4s ",TransX[i]);
142: for (i=MM; i<M; i++) fprintf(fp," %4s ",TransX[i]);
143: for (i=NN; i<N; i++) fprintf(fp," %4s ",TransX[i]);
144: for (i=CC; i<C; i++) fprintf(fp," %4s ",TransD[i]);
145: for (i=LL; i<L; i++) fprintf(fp," %4s ",TransD[i]);
146: for (i=MM; i<M; i++) fprintf(fp," %4s ",TransD[i]);
147: for (i=NN; i<N; i++) fprintf(fp," %4s ",TransD[i]);
148: fprintf(fp,"\n");
149: }else {
150: fprintf(fp,"Unknown graduation variable specification.\n\n");
151: }
152: fprintf(fp,"The homogenization variable is : ");
153: fprintf(fp," %4s ",TransD[0]);
154: fprintf(fp,"\n");
155:
156:
157:
158: fprintf(fp,"-------------------------------------------\n");
159: fprintf(fp,"Output order : ");
160: for (i=0; i<2*N; i++) {
161: if (ringp->outputOrder[i] < N) {
162: fprintf(fp,"%s ",TransX[ringp->outputOrder[i]]);
163: }else{
164: fprintf(fp,"%s ",TransD[(ringp->outputOrder[i])-N]);
165: }
166: }
167: fprintf(fp,"\n");
168:
169: if (ringp->multiplication == mpMult_poly) {
170: mtype = "poly";
171: }else if (ringp->multiplication == mpMult_diff) {
172: mtype = "diff";
173: }else if (ringp->multiplication == mpMult_difference) {
174: mtype = "difference";
175: }else {
176: mtype = "unknown";
177: }
1.16 takayama 178: fprintf(fp,"Multiplication function --%s(%p).\n",
179: mtype, ringp->multiplication);
1.1 maekawa 180: if (ringp->schreyer) {
181: fprintf(fp,"schreyer=1, gbListTower=");
182: printObjectList((struct object *)(ringp->gbListTower));
183: fprintf(fp,"\n");
184: }
1.7 takayama 185: if (ringp->degreeShiftSize) {
1.8 takayama 186: fprintf(fp,"degreeShift vector (N=%d,Size=%d)= \n[\n",ringp->degreeShiftN,ringp->degreeShiftSize);
1.7 takayama 187: {
1.8 takayama 188: int i,j;
189: for (i=0; i<ringp->degreeShiftN; i++) {
190: fprintf(fp," [");
191: for (j=0; j< ringp->degreeShiftSize; j++) {
192: fprintf(fp," %d ",ringp->degreeShift[i*(ringp->degreeShiftSize)+j]);
193: }
194: fprintf(fp,"]\n");
1.7 takayama 195: }
196: }
197: fprintf(fp,"]\n");
198: }
199: fprintf(fp,"--- weight vectors ---\n");
1.1 maekawa 200: if (level) printOrder(ringp);
1.13 takayama 201:
202: if (ringp->partialEcart) {
203: fprintf(fp,"--- partialEcartGlobalVarX ---\n");
204: for (i=0; i<ringp->partialEcart; i++) {
205: fprintf(fp," %4s ",TransX[ringp->partialEcartGlobalVarX[i]]);
206: }
207: fprintf(fp,"\n");
208: }
1.1 maekawa 209:
210: if (ringp->next != (struct ring *)NULL) {
211: fprintf(fp,"\n\n-------- The next ring is .... --------------\n");
212: showRing(level,ringp->next);
213: }
214: }
215:
216: /***************************************************************
217: functions related to order
218: ******************************************************************/
219: #define xtoi(k) ((N-1)-(k))
220: #define dtoi(k) ((2*N-1)-(k))
221: #define itox(k) ((N-1)-(k))
222: #define itod(k) ((2*N-1)-(k))
223: #define isX(i) (i<N? 1: 0)
224: #define isD(i) (i<N? 0: 1)
225: /****************************************************
226: i : 0 1 N-1 N 2N-1
227: x :x_{N-1} x_{N-2} x_0
228: d : D_{N-1} D_{0}
229: if (isX(i)) x_{itox(i)}
230: if (isD(i)) D_{itod(i)}
231: ******************************************************/
232: /* xtoi(0):N-1 xtoi(1):N-2 ....
233: dtoi(0):2N-1 dtoi(1):2N-2 ...
234: itod(N):N-1 dtoi(N-1):N ...
235: */
236:
237: void printOrder(ringp)
1.4 takayama 238: struct ring *ringp;
1.1 maekawa 239: {
240: int i,j;
241: FILE *fp;
242: char tmp[100];
243: int N,M,L,C,NN,MM,LL,CC;
244: char **TransX,**TransD;
245: int *Order;
246: int P;
247: int omsize;
248: extern char *F_isSameComponent;
249:
250: N=ringp->n; M = ringp->m; L = ringp->l; C = ringp->c;
251: NN=ringp->nn; MM = ringp->mm; LL = ringp->ll; CC = ringp->cc;
252: TransX = ringp->x; TransD = ringp->D;
253: Order = ringp->order;
254: P = ringp->p;
255: omsize = ringp->orderMatrixSize;
256:
257: fp = stdout;
258:
259:
260: for (i=0; i<2*N; i++) printf("%4d",i);
261: fprintf(fp,"\n");
262:
263: /* print variables names */
264: for (i=0; i<N; i++) {
265: sprintf(tmp,"x%d",N-1-i);
266: fprintf(fp,"%4s",tmp);
267: }
268: for (i=0; i<N; i++) {
269: sprintf(tmp,"D%d",N-1-i);
270: fprintf(fp,"%4s",tmp);
271: }
272: fprintf(fp,"\n");
273:
274: /* print identifier names */
275: for (i=0; i<N; i++) fprintf(fp,"%4s",TransX[itox(i)]);
276: for (i=N; i<2*N; i++) fprintf(fp,"%4s",TransD[itod(i)]);
277: fprintf(fp,"\n");
278:
279: /* print D: differential DE: differential, should be eliminated
1.4 takayama 280: E: difference
281: Q: q-difference
282: C: commutative
1.1 maekawa 283: */
284: if (strcmp(F_isSameComponent,"x")== 0 || strcmp(F_isSameComponent,"xd")==0) {
285: for (i=0; i<N; i++) {
286: if ((NN<=itox(i)) && (itox(i)<N)) fprintf(fp,"%4s","DE");
287: if ((M<=itox(i)) && (itox(i)<NN)) fprintf(fp,"%4s","D");
288: if ((MM<=itox(i)) && (itox(i)<M)) fprintf(fp,"%4s","EE");
289: if ((L<=itox(i)) && (itox(i)<MM)) fprintf(fp,"%4s","E");
290: if ((LL<=itox(i)) && (itox(i)<L)) fprintf(fp,"%4s","QE");
291: if ((C<=itox(i)) && (itox(i)<LL)) fprintf(fp,"%4s","Q");
292: if ((CC<=itox(i)) && (itox(i)<C)) fprintf(fp,"%4s","CE");
293: if ((0<=itox(i)) && (itox(i)<CC)) fprintf(fp,"%4s","C");
294: }
295: }
296: if (strcmp(F_isSameComponent,"x")==0) {
297: for (i=N; i<2*N; i++) {
298: if ((M<=itod(i)) && (itod(i)<N)) fprintf(fp,"%4s","D");
299: if ((L<=itod(i)) && (itod(i)<M)) fprintf(fp,"%4s","E");
300: if ((C<=itod(i)) && (itod(i)<L)) fprintf(fp,"%4s","Q");
301: if ((0<=itod(i)) && (itod(i)<C)) fprintf(fp,"%4s","C");
302: }
303: }else if (strcmp(F_isSameComponent,"xd")==0) {
304: for (i=N; i<2*N; i++) {
305: if ((NN<=itod(i)) && (itod(i)<N)) fprintf(fp,"%4s","DE");
306: if ((M<=itod(i)) && (itod(i)<NN)) fprintf(fp,"%4s","D");
307: if ((MM<=itod(i)) && (itod(i)<M)) fprintf(fp,"%4s","EE");
308: if ((L<=itod(i)) && (itod(i)<MM)) fprintf(fp,"%4s","E");
309: if ((LL<=itod(i)) && (itod(i)<L)) fprintf(fp,"%4s","QE");
310: if ((C<=itod(i)) && (itod(i)<LL)) fprintf(fp,"%4s","Q");
311: if ((CC<=itod(i)) && (itod(i)<C)) fprintf(fp,"%4s","CE");
312: if ((0<=itod(i)) && (itod(i)<CC)) fprintf(fp,"%4s","C");
313: }
314: } else {
315: fprintf(fp,"Unknown graduation variable type.\n");
316: }
317: fprintf(fp,"\n");
318:
319: for (i=0; i< omsize; i++) {
320: for (j=0; j<2*N; j++) {
321: fprintf(fp,"%4d", Order[i*2*N+j]);
322: }
323: fprintf(fp,"\n");
324: }
325: fprintf(fp,"\n");
326:
327: }
328:
329: struct object oGetOrderMatrix(struct ring *ringp)
330: {
1.14 takayama 331: struct object rob = OINIT;
332: struct object ob2 = OINIT;
1.1 maekawa 333: int n,i,j,m;
334: int *om;
335: n = ringp->n;
336: m = ringp->orderMatrixSize;
337: om = ringp->order;
338: if (m<=0) m = 1;
339: rob = newObjectArray(m);
340: for (i=0; i<m; i++) {
341: ob2 = newObjectArray(2*n);
342: for (j=0; j<2*n; j++) {
343: putoa(ob2,j,KpoInteger(om[2*n*i+j]));
344: }
345: putoa(rob,i,ob2);
346: }
347: return(rob);
348: }
349:
350:
351: int mmLarger_matrix(ff,gg)
1.4 takayama 352: POLY ff; POLY gg;
1.1 maekawa 353: {
354: int exp[2*N0]; /* exponents */
355: int i,k;
356: int sum,flag;
357: int *Order;
358: int N;
359: MONOMIAL f,g;
360: struct ring *rp;
361: int in2;
362: int *from, *to;
363: int omsize;
1.7 takayama 364: int dssize;
1.8 takayama 365: int dsn;
1.7 takayama 366: int *degreeShiftVector;
1.1 maekawa 367:
368: if (ff == POLYNULL ) {
369: if (gg == POLYNULL) return( 2 );
370: else return( 0 );
371: }
372: if (gg == POLYNULL) {
373: if (ff == POLYNULL) return( 2 );
374: else return( 1 );
375: }
376: f = ff->m; g=gg->m;
377:
378: rp = f->ringp;
379: Order = rp->order;
380: N = rp->n;
381: from = rp->from;
382: to = rp->to;
383: omsize = rp->orderMatrixSize;
1.7 takayama 384: if (dssize = rp->degreeShiftSize) {
385: degreeShiftVector = rp->degreeShift; /* Note. 2003.06.26 */
1.8 takayama 386: dsn = rp->degreeShiftN;
1.7 takayama 387: }
1.1 maekawa 388:
389: flag = 1;
390: for (i=N-1,k=0; i>=0; i--,k++) {
391: exp[k] = (f->e[i].x) - (g->e[i].x);
392: exp[k+N] = (f->e[i].D) - (g->e[i].D);
393: if ((exp[k] != 0) || (exp[k+N] != 0)) flag =0;
394: }
395: if (flag==1) return(2);
396: /* exp > 0 <---> f>g
397: exp = 0 <---> f=g
398: exp < 0 <---> f<g
399: */
400: for (i=0; i< omsize; i++) {
401: sum = 0; in2 = i*2*N;
402: /* for (k=0; k<2*N; k++) sum += exp[k]*Order[in2+k]; */
403: for (k=from[i]; k<to[i]; k++) sum += exp[k]*Order[in2+k];
1.8 takayama 404: if (dssize && ( i < dsn)) { /* Note, 2003.06.26 */
1.7 takayama 405: if ((f->e[N-1].x < dssize) && (f->e[N-1].x >= 0) &&
406: (g->e[N-1].x < dssize) && (g->e[N-1].x >= 0)) {
1.8 takayama 407: sum += degreeShiftVector[i*dssize+ (f->e[N-1].x)]
408: -degreeShiftVector[i*dssize+ (g->e[N-1].x)];
1.7 takayama 409: }else{
1.9 takayama 410: /*warningOrder("Size mismatch in the degree shift vector. It is ignored.");*/
1.7 takayama 411: }
412: }
1.1 maekawa 413: if (sum > 0) return(1);
414: if (sum < 0) return(0);
415: }
416: return(2);
417: }
418:
419: /* This should be used in case of q */
420: int mmLarger_qmatrix(ff,gg)
1.4 takayama 421: POLY ff; POLY gg;
1.1 maekawa 422: {
423: int exp[2*N0]; /* exponents */
424: int i,k;
425: int sum,flag;
426: int *Order;
427: int N;
428: MONOMIAL f,g;
429: int omsize;
430:
431: if (ff == POLYNULL ) {
432: if (gg == POLYNULL) return( 2 );
433: else return( 0 );
434: }
435: if (gg == POLYNULL) {
436: if (ff == POLYNULL) return( 2 );
437: else return( 1 );
438: }
439: f = ff->m; g = gg->m;
440: Order = f->ringp->order;
441: N = f->ringp->n;
442: omsize = f->ringp->orderMatrixSize;
443:
444: flag = 1;
445: for (i=N-1,k=0; i>=0; i--,k++) {
446: exp[k] = (f->e[i].x) - (g->e[i].x);
447: exp[k+N] = (f->e[i].D) - (g->e[i].D);
448: if ((exp[k] != 0) || (exp[k+N] != 0)) flag =0;
449: }
450: if (flag==1) return(2);
451: /* exp > 0 <---> f>g
452: exp = 0 <---> f=g
453: exp < 0 <---> f<g
454: */
455: for (i=0; i< omsize; i++) {
456: sum = 0;
457: /* In case of q, you should do as follows */
458: for (k=0; k<N-1; k++) sum += exp[k]*Order[i*2*N+k]; /* skip k= N-1 -->q */
459: for (k=N; k<2*N-1; k++) sum += exp[k]*Order[i*2*N+k]; /* SKip k= 2*N-1 */
460: if (sum > 0) return(1);
461: else if (sum < 0) return(0);
462: }
463: if (exp[N-1] > 0) return(1);
464: else if (exp[N-1] < 0) return(0);
465: else return(2);
466: }
467:
468: /* x(N-1)>x(N-2)>....>D(N-1)>....>D(0) */
1.17 ! takayama 469: int mmLarger_pureLexicographic(f,g)
1.4 takayama 470: POLY f;
471: POLY g;
1.1 maekawa 472: {
473: int i,r;
474: int n;
475: MONOMIAL fm,gm;
476: /* Note that this function ignores the order matrix of the given
477: ring. */
478: if (f == POLYNULL ) {
479: if (g == POLYNULL) return( 2 );
480: else return( 0 );
481: }
482: if (g == POLYNULL) {
483: if (f == POLYNULL) return( 2 );
484: else return( 1 );
485: }
486:
487:
488: fm = f->m; gm = g->m;
489: n = fm->ringp->n;
490: for (i=n-1; i>=0; i--) {
491: r = (fm->e[i].x) - (gm->e[i].x);
492: if (r > 0) return(1);
493: else if (r < 0) return(0);
494: else ;
495: }
496:
497: for (i=n-1; i>=0; i--) {
498: r = (fm->e[i].D) - (gm->e[i].D);
499: if (r > 0) return(1);
500: else if (r < 0) return(0);
501: else ;
502: }
503:
504: return(2);
505:
506: }
507:
508:
509: void setFromTo(ringp)
1.4 takayama 510: struct ring *ringp;
1.1 maekawa 511: {
512: int n;
513: int i,j,oasize;
514: if (ringp->order == (int *)NULL) errorOrder("setFromTo(); no order matrix.");
515: n = (ringp->n)*2;
516: oasize = ringp->orderMatrixSize;
517: ringp->from = (int *)sGC_malloc(sizeof(int)*oasize);
518: ringp->to = (int *)sGC_malloc(sizeof(int)*oasize);
519: if (ringp->from == (int *)NULL || ringp->to == (int *)NULL) {
520: errorOrder("setFromTo(): No memory.");
521: }
522: for (i=0; i<oasize; i++) {
523: ringp->from[i] = 0; ringp->to[i] = n;
524: for (j=0; j<n; j++) {
525: if (ringp->order[i*n+j] != 0) {
1.4 takayama 526: ringp->from[i] = j;
527: break;
1.1 maekawa 528: }
529: }
530: for (j=n-1; j>=0; j--) {
531: if (ringp->order[i*n+j] != 0) {
1.4 takayama 532: ringp->to[i] = j+1;
533: break;
1.1 maekawa 534: }
535: }
536: }
537: }
538:
539: /* It ignores h and should be used with mmLarger_tower */
540: /* cf. mmLarger_matrix. h always must be checked at last. */
541: static int mmLarger_matrix_schreyer(ff,gg)
1.4 takayama 542: POLY ff; POLY gg;
1.1 maekawa 543: {
544: int exp[2*N0]; /* exponents */
545: int i,k;
546: int sum,flag;
547: int *Order;
548: int N;
549: MONOMIAL f,g;
550: struct ring *rp;
551: int in2;
552: int *from, *to;
553: int omsize;
554:
555: if (ff == POLYNULL ) {
556: if (gg == POLYNULL) return( 2 );
557: else return( 0 );
558: }
559: if (gg == POLYNULL) {
560: if (ff == POLYNULL) return( 2 );
561: else return( 1 );
562: }
563: f = ff->m; g=gg->m;
564:
565: rp = f->ringp;
566: Order = rp->order;
567: N = rp->n;
568: from = rp->from;
569: to = rp->to;
570: omsize = rp->orderMatrixSize;
571:
572: flag = 1;
573: for (i=N-1,k=0; i>0; i--,k++) {
574: exp[k] = (f->e[i].x) - (g->e[i].x);
575: exp[k+N] = (f->e[i].D) - (g->e[i].D);
576: if ((exp[k] != 0) || (exp[k+N] != 0)) flag =0;
577: }
578: exp[N-1] = (f->e[0].x) - (g->e[0].x);
579: exp[2*N-1] = 0; /* f->e[0].D - g->e[0].D. Ignore h! */
580: if ((exp[N-1] != 0) || (exp[2*N-1] != 0)) flag =0;
581:
582: if (flag==1) return(2);
583: /* exp > 0 <---> f>g
584: exp = 0 <---> f=g
585: exp < 0 <---> f<g
586: */
587: for (i=0; i< omsize; i++) {
588: sum = 0; in2 = i*2*N;
589: /* for (k=0; k<2*N; k++) sum += exp[k]*Order[in2+k]; */
590: for (k=from[i]; k<to[i]; k++) sum += exp[k]*Order[in2+k];
591: if (sum > 0) return(1);
592: if (sum < 0) return(0);
593: }
594: return(2);
595: }
596:
597: int mmLarger_tower(POLY f,POLY g) {
598: struct object *gbList;
599: int r;
600: if (f == POLYNULL) {
601: if (g == POLYNULL) return(2);
602: else return(0);
603: }
604: if (g == POLYNULL) {
605: if (f == POLYNULL) return(2);
606: else return(1);
607: }
608: if (!(f->m->ringp->schreyer) || !(g->m->ringp->schreyer))
609: return(mmLarger_matrix(f,g));
1.4 takayama 610: /* modifiable: mmLarger_qmatrix */
1.1 maekawa 611: gbList = (struct object *)(g->m->ringp->gbListTower);
612: if (gbList == NULL) return(mmLarger_matrix(f,g));
1.4 takayama 613: /* modifiable: mmLarger_qmatrix */
1.1 maekawa 614: if (gbList->tag != Slist) {
615: warningOrder("mmLarger_tower(): gbList must be in Slist.\n");
616: return(1);
617: }
618: if (klength(gbList) ==0) return(mmLarger_matrix(f,g));
1.4 takayama 619: /* modifiable: mmLarger_qmatrix */
1.1 maekawa 620:
621: r = mmLarger_tower3(f,g,gbList);
622: /* printf("mmLarger_tower3(%s,%s) --> %d\n",POLYToString(head(f),'*',1),POLYToString(head(g),'*',1),r); */
623: if (r == 2) { /* Now, compare by h */
624: if (f->m->e[0].D > g->m->e[0].D) return(1);
625: else if (f->m->e[0].D < g->m->e[0].D) return(0);
626: else return(2);
627: }else{
628: return(r);
629: }
630: }
631:
632: int mmLarger_tower3(POLY f,POLY g,struct object *gbList)
633: { /* gbList is assumed to be Slist */
634: int n,fv,gv,t,r,nn;
635: POLY fm;
636: POLY gm;
1.14 takayama 637: struct object gb = OINIT;
1.1 maekawa 638:
639: if (f == POLYNULL) {
640: if (g == POLYNULL) return(2);
641: else return(0);
642: }
643: if (g == POLYNULL) {
644: if (f == POLYNULL) return(2);
645: else return(1); /* It assumes the zero is the minimum element!! */
646: }
647: n = f->m->ringp->n;
648: nn = f->m->ringp->nn;
649: /* critical and modifiable */ /* m e_u > m e_v <==> m g_u > m g_v */
1.4 takayama 650: /* or equal and u < v */
1.1 maekawa 651: fv = f->m->e[nn].x ; /* extract component (vector) number of f! */
652: gv = g->m->e[nn].x ;
653: if (fv == gv) { /* They have the same component number. */
654: return(mmLarger_matrix_schreyer(f,g));
655: }
656:
657: if (gbList == NULL) return(mmLarger_matrix_schreyer(f,g));
1.4 takayama 658: /* modifiable: mmLarger_qmatrix */
1.1 maekawa 659: if (gbList->tag != Slist) {
660: warningOrder("mmLarger_tower(): gbList must be in Slist.\n");
661: return(1);
662: }
663: if (klength(gbList) ==0) return(mmLarger_matrix(f,g));
1.4 takayama 664: /* modifiable: mmLarger_qmatrix */
1.1 maekawa 665: gb = car(gbList); /* each entry must be monomials */
666: if (gb.tag != Sarray) {
667: warningOrder("mmLarger_tower3(): car(gbList) must be an array.\n");
668: return(1);
669: }
670: t = getoaSize(gb);
671: if (t == 0) return(mmLarger_tower3(f,g,cdr(gbList)));
672:
673: fm = pmCopy(head(f)); fm->m->e[nn].x = 0; /* f is not modified. */
674: gm = pmCopy(head(g)); gm->m->e[nn].x = 0;
675: if (fv >= t || gv >= t) {
676: warningOrder("mmLarger_tower3(): incompatible input and gbList.\n");
677: printf("Length of gb is %d, f is %s, g is %s\n",t,KPOLYToString(f),
1.4 takayama 678: KPOLYToString(g));
1.3 takayama 679: KSexecuteString(" show_ring ");
1.1 maekawa 680: return(1);
681: }
682: /* mpMult_poly is too expensive to call. @@@*/
683: r = mmLarger_tower3(mpMult_poly(fm,KopPOLY(getoa(gb,fv))),
684: mpMult_poly(gm,KopPOLY(getoa(gb,gv))),
685: cdr(gbList));
686: if (r != 2) return(r);
687: else if (fv == gv) return(2);
688: else if (fv > gv) return(0); /* modifiable */
689: else if (fv < gv) return(1); /* modifiable */
690: }
1.11 takayama 691:
692: static struct object auxPruneZeroRow(struct object ob) {
693: int i,m,size;
1.14 takayama 694: struct object obt = OINIT;
695: struct object rob = OINIT;
1.11 takayama 696: m = getoaSize(ob);
697: size=0;
698: for (i=0; i<m; i++) {
699: obt = getoa(ob,i);
700: if (getoaSize(obt) != 0) size++;
701: }
702: if (size == m) return ob;
703: rob = newObjectArray(size);
704: for (i=0, size=0; i<m; i++) {
705: obt = getoa(ob,i);
706: if (getoaSize(obt) != 0) {
707: putoa(rob,size,obt); size++;
708: }
709: }
710: return rob;
711: }
1.12 takayama 712: static struct object oRingToOXringStructure_long(struct ring *ringp)
1.10 takayama 713: {
1.14 takayama 714: struct object rob = OINIT;
715: struct object ob2 = OINIT;
716: struct object obMat = OINIT;
717: struct object obV = OINIT;
718: struct object obShift = OINIT;
719: struct object obt = OINIT;
1.10 takayama 720: char **TransX; char **TransD;
721: int n,i,j,m,p,nonzero;
722: int *om;
723: n = ringp->n;
724: m = ringp->orderMatrixSize;
725: om = ringp->order;
726: TransX = ringp->x; TransD = ringp->D;
727: if (m<=0) m = 1;
728: /*test: (1). getRing /rr set rr (oxRingStructure) dc */
729: obMat = newObjectArray(m);
730: for (i=0; i<m; i++) {
731: nonzero = 0;
732: for (j=0; j<2*n; j++) {
733: if (om[2*n*i+j] != 0) nonzero++;
734: }
735: ob2 = newObjectArray(nonzero*2);
736: nonzero=0;
737: for (j=0; j<2*n; j++) {
738: /* fprintf(stderr,"%d, ",nonzero); */
739: if (om[2*n*i+j] != 0) {
740: if (j < n) {
741: putoa(ob2,nonzero,KpoString(TransX[n-1-j])); nonzero++;
742: }else{
743: putoa(ob2,nonzero,KpoString(TransD[n-1-(j-n)])); nonzero++;
744: }
745: putoa(ob2,nonzero,KpoUniversalNumber(newUniversalNumber(om[2*n*i+j]))); nonzero++;
746: }
747: }
748: /* printObject(ob2,0,stderr); fprintf(stderr,".\n"); */
749: putoa(obMat,i,ob2);
750: }
1.11 takayama 751: obMat = auxPruneZeroRow(obMat);
1.10 takayama 752: /* printObject(obMat,0,stderr); */
753:
754: obV = newObjectArray(2*n);
755: for (i=0; i<n; i++) putoa(obV,i,KpoString(TransX[n-1-i]));
756: for (i=0; i<n; i++) putoa(obV,i+n,KpoString(TransD[n-1-i]));
757: /* printObject(obV,0,stderr); */
758:
759: if (ringp->degreeShiftSize) {
760: /*test:
761: [(x) ring_of_differential_operators [[(x)]] weight_vector 0
762: [(weightedHomogenization) 1 (degreeShift) [[1 2 1]]] ] define_ring ;
763: (1). getRing /rr set rr (oxRingStructure) dc message
764: */
765: obShift = newObjectArray(ringp->degreeShiftN);
766: for (i=0; i<ringp->degreeShiftN; i++) {
767: obt = newObjectArray(ringp->degreeShiftSize);
768: for (j=0; j< ringp->degreeShiftSize; j++) {
769: putoa(obt,j,KpoUniversalNumber(newUniversalNumber(ringp->degreeShift[i*(ringp->degreeShiftSize)+j])));
770: }
771: putoa(obShift,i,obt);
772: }
773: /* printObject(obShift,0,stderr); */
774: }
775:
776: p = 0;
777: if (ringp->degreeShiftSize) {
778: rob = newObjectArray(3);
779: obt = newObjectArray(2);
780: putoa(obt,0,KpoString("degreeShift"));
781: putoa(obt,1,obShift);
782: putoa(rob,p, obt); p++;
783: }else {
784: rob = newObjectArray(2);
785: }
786:
787: obt = newObjectArray(2);
788: putoa(obt,0,KpoString("v"));
789: putoa(obt,1,obV);
790: putoa(rob,p, obt); p++;
791:
792: obt = newObjectArray(2);
793: putoa(obt,0,KpoString("order"));
794: putoa(obt,1,obMat);
795: putoa(rob,p, obt); p++;
796:
1.12 takayama 797: return(rob);
798: }
799: static int auxEffectiveVar(int idx,int n) {
800: int x;
801: if (idx < n) x=1; else x=0;
802: if (x) {
803: if ((idx >= 1) && (idx < n-1)) return 1;
804: else return 0;
805: }else{
806: if ( 1 <= idx-n ) return 1;
807: else return 0;
808: }
809: }
810: /*test:
811: [(x,y) ring_of_differential_operators [[(Dx) 1 (Dy) 1]]
812: weight_vector 0] define_ring
813: (x). getRing (oxRingStructure) dc ::
814: */
815: static struct object oRingToOXringStructure_short(struct ring *ringp)
816: {
1.14 takayama 817: struct object rob = OINIT;
818: struct object ob2 = OINIT;
819: struct object obMat = OINIT;
820: struct object obV = OINIT;
821: struct object obShift = OINIT;
822: struct object obt = OINIT;
1.12 takayama 823: char **TransX; char **TransD;
824: int n,i,j,m,p,nonzero;
825: int *om;
826: n = ringp->n;
827: m = ringp->orderMatrixSize;
828: om = ringp->order;
829: TransX = ringp->x; TransD = ringp->D;
830: if (m<=0) m = 1;
831: /*test: (1). getRing /rr set rr (oxRingStructure) dc */
832: obMat = newObjectArray(m);
833: for (i=0; i<m; i++) {
834: nonzero = 0;
835: for (j=0; j<2*n; j++) {
836: if ((om[2*n*i+j] != 0) && auxEffectiveVar(j,n)) nonzero++;
837: }
838: ob2 = newObjectArray(nonzero*2);
839: nonzero=0;
840: for (j=0; j<2*n; j++) {
841: /* fprintf(stderr,"%d, ",nonzero); */
842: if ((om[2*n*i+j] != 0) && auxEffectiveVar(j,n)) {
843: if (j < n) {
844: putoa(ob2,nonzero,KpoString(TransX[n-1-j])); nonzero++;
845: }else{
846: putoa(ob2,nonzero,KpoString(TransD[n-1-(j-n)])); nonzero++;
847: }
848: putoa(ob2,nonzero,KpoUniversalNumber(newUniversalNumber(om[2*n*i+j]))); nonzero++;
849: }
850: }
851: /* printObject(ob2,0,stderr); fprintf(stderr,".\n"); */
852: putoa(obMat,i,ob2);
853: }
854: obMat = auxPruneZeroRow(obMat);
855: /* printObject(obMat,0,stderr); */
856:
857: obV = newObjectArray(2*n-3);
858: for (i=0; i<n-2; i++) putoa(obV,i,KpoString(TransX[n-1-i-1]));
859: for (i=0; i<n-1; i++) putoa(obV,i+n-2,KpoString(TransD[n-1-i-1]));
860: /* printObject(obV,0,stderr); */
861:
862: if (ringp->degreeShiftSize) {
863: /*test:
864: [(x) ring_of_differential_operators [[(x)]] weight_vector 0
865: [(weightedHomogenization) 1 (degreeShift) [[1 2 1]]] ] define_ring ;
866: (1). getRing /rr set rr (oxRingStructure) dc message
867: */
868: obShift = newObjectArray(ringp->degreeShiftN);
869: for (i=0; i<ringp->degreeShiftN; i++) {
870: obt = newObjectArray(ringp->degreeShiftSize);
871: for (j=0; j< ringp->degreeShiftSize; j++) {
872: putoa(obt,j,KpoUniversalNumber(newUniversalNumber(ringp->degreeShift[i*(ringp->degreeShiftSize)+j])));
873: }
874: putoa(obShift,i,obt);
875: }
876: /* printObject(obShift,0,stderr); */
877: }
878:
879: p = 0;
880: if (ringp->degreeShiftSize) {
881: rob = newObjectArray(3);
882: obt = newObjectArray(2);
883: putoa(obt,0,KpoString("degreeShift"));
884: putoa(obt,1,obShift);
885: putoa(rob,p, obt); p++;
886: }else {
887: rob = newObjectArray(2);
888: }
889:
890: obt = newObjectArray(2);
891: putoa(obt,0,KpoString("v"));
892: putoa(obt,1,obV);
893: putoa(rob,p, obt); p++;
894:
895: obt = newObjectArray(2);
896: putoa(obt,0,KpoString("order"));
897: putoa(obt,1,obMat);
898: putoa(rob,p, obt); p++;
899:
900: return(rob);
901: }
902: struct object oRingToOXringStructure(struct ring *ringp)
903: {
1.14 takayama 904: struct object rob = OINIT;
905: struct object tob = OINIT;
1.12 takayama 906: rob = newObjectArray(2);
907: tob = oRingToOXringStructure_short(ringp);
908: putoa(rob,0,tob);
909: tob = oRingToOXringStructure_long(ringp);
910: putoa(rob,1,tob);
1.10 takayama 911: return(rob);
912: }
913:
1.1 maekawa 914: static void warningOrder(s)
1.4 takayama 915: char *s;
1.1 maekawa 916: {
917: fprintf(stderr,"Warning in order.c: %s\n",s);
918: }
919:
920: static void errorOrder(s)
1.4 takayama 921: char *s;
1.1 maekawa 922: {
923: fprintf(stderr,"order.c: %s\n",s);
924: exit(14);
925: }
926:
927:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>