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