![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to C.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / engine|
Return to C.c CVS log | Up to [local] / OpenXM_contrib2 / asir2000 / engine |
1.1 noro 1: /* $OpenXM: OpenXM/src/asir99/engine/C.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
2: #include "ca.h"
3: #include "inline.h"
4: #include "base.h"
5:
6: V up_var;
7:
8: /* binary has at least 32 leading 0 chars. */
9: void binaryton(binary,np)
10: char *binary;
11: N *np;
12: {
13: int i,w,len;
14: N n;
15: char buf[33];
16:
17: binary += strlen(binary)%32;
18: len = strlen(binary);
19: w = len/32; /* sufficient for holding binary */
20: n = NALLOC(w);
21: for ( i = 0; i < w; i++ ) {
22: strncpy(buf,binary+len-32*(i+1),32); buf[32] = 0;
23: n->b[i] = strtoul(buf,0,2);
24: }
25: for ( i = w-1; i >= 0 && !n->b[i]; i-- );
26: if ( i < 0 )
27: *np = 0;
28: else {
29: n->p = i+1;
30: *np = n;
31: }
32: }
33:
34: /* hex has at least 8 leading 0 chars. */
35: void hexton(hex,np)
36: char *hex;
37: N *np;
38: {
39: int i,w,len;
40: N n;
41: char buf[9];
42:
43: hex += strlen(hex)%8;
44: len = strlen(hex);
45: w = len/8; /* sufficient for holding hex */
46: n = NALLOC(w);
47: for ( i = 0; i < w; i++ ) {
48: strncpy(buf,hex+len-8*(i+1),8); buf[8] = 0;
49: n->b[i] = strtoul(buf,0,16);
50: }
51: for ( i = w-1; i >= 0 && !n->b[i]; i-- );
52: if ( i < 0 )
53: *np = 0;
54: else {
55: n->p = i+1;
56: *np = n;
57: }
58: }
59:
60: void ntobn(base,n,nrp)
61: int base;
62: N n,*nrp;
63: {
64: int i,d,plc;
65: unsigned int *c,*x,*w;
66: unsigned int r;
67: L m;
68: N nr;
69:
70: if ( !n ) {
71: *nrp = NULL;
72: return;
73: }
74:
75: d = PL(n);
76: w = BD(n);
77:
78: for ( i = 1, m = 1; m <= LBASE/(L)base; m *= base, i++ );
79:
80: c = (unsigned int *)W_ALLOC(d*i+1);
81: x = (unsigned int *)W_ALLOC(d+1);
82: for ( i = 0; i < d; i++ )
83: x[i] = w[i];
84: for ( plc = 0; d >= 1; plc++ ) {
85: for ( i = d - 1, r = 0; i >= 0; i-- ) {
86: DSAB((unsigned int)base,r,x[i],x[i],r)
87: }
88: c[plc] = r;
89: if ( !x[d-1] ) d--;
90: }
91:
92: *nrp = nr = NALLOC(plc); INITRC(nr);
93: PL(nr) = plc;
94: for ( i = 0; i < plc; i++ )
95: BD(nr)[i] = c[i];
96: }
97:
98: void bnton(base,n,nrp)
99: int base;
100: N n,*nrp;
101: {
102: unsigned int carry;
103: unsigned int *x,*w;
104: int i,j,d,plc;
105: N nr;
106:
107: if ( !n ) {
108: *nrp = 0;
109: return;
110: }
111:
112: d = PL(n);
113: w = BD(n);
114: x = (unsigned int *)W_ALLOC(d + 1);
115:
116: for ( plc = 0, i = d - 1; i >= 0; i-- ) {
117: for ( carry = w[i],j = 0; j < plc; j++ ) {
118: DMA(x[j],(unsigned int)base,carry,carry,x[j])
119: }
120: if ( carry ) x[plc++] = carry;
121: }
122: *nrp = nr = NALLOC(plc); INITRC(nr);
123: PL(nr) = plc;
124: for ( i = 0; i < plc; i++ )
125: BD(nr)[i] = x[i];
126: }
127:
128: void ptomp(m,p,pr)
129: int m;
130: P p;
131: P *pr;
132: {
133: DCP dc,dcr,dcr0;
134: Q q;
135: unsigned int a,b;
136: P t;
137: MQ s;
138:
139: if ( !p )
140: *pr = 0;
141: else if ( NUM(p) ) {
142: q = (Q)p;
143: a = rem(NM(q),m);
144: if ( a && (SGN(q) < 0) )
145: a = m-a;
146: b = !DN(q)?1:rem(DN(q),m);
147: if ( !b )
148: error("ptomp : denominator = 0");
149: a = dmar(a,invm(b,m),0,m); STOMQ(a,s); *pr = (P)s;
150: } else {
151: for ( dc = DC(p), dcr0 = 0; dc; dc = NEXT(dc) ) {
152: ptomp(m,COEF(dc),&t);
153: if ( t ) {
154: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = t;
155: }
156: }
157: if ( !dcr0 )
158: *pr = 0;
159: else {
160: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
161: }
162: }
163: }
164:
165: void mptop(f,gp)
166: P f;
167: P *gp;
168: {
169: DCP dc,dcr,dcr0;
170: Q q;
171:
172: if ( !f )
173: *gp = 0;
174: else if ( NUM(f) )
175: STOQ(CONT((MQ)f),q),*gp = (P)q;
176: else {
177: for ( dc = DC(f), dcr0 = 0; dc; dc = NEXT(dc) ) {
178: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); mptop(COEF(dc),&COEF(dcr));
179: }
180: NEXT(dcr) = 0; MKP(VR(f),dcr0,*gp);
181: }
182: }
183:
184: void ptolmp(p,pr)
185: P p;
186: P *pr;
187: {
188: DCP dc,dcr,dcr0;
189: LM a;
190: P t;
191:
192: if ( !p )
193: *pr = 0;
194: else if ( NUM(p) ) {
195: qtolm((Q)p,&a); *pr = (P)a;
196: } else {
197: for ( dc = DC(p), dcr0 = 0; dc; dc = NEXT(dc) ) {
198: ptolmp(COEF(dc),&t);
199: if ( t ) {
200: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); COEF(dcr) = t;
201: }
202: }
203: if ( !dcr0 )
204: *pr = 0;
205: else {
206: NEXT(dcr) = 0; MKP(VR(p),dcr0,*pr);
207: }
208: }
209: }
210:
211: void lmptop(f,gp)
212: P f;
213: P *gp;
214: {
215: DCP dc,dcr,dcr0;
216: Q q;
217:
218: if ( !f )
219: *gp = 0;
220: else if ( NUM(f) ) {
221: NTOQ(((LM)f)->body,1,q); *gp = (P)q;
222: } else {
223: for ( dc = DC(f), dcr0 = 0; dc; dc = NEXT(dc) ) {
224: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc); lmptop(COEF(dc),&COEF(dcr));
225: }
226: NEXT(dcr) = 0; MKP(VR(f),dcr0,*gp);
227: }
228: }
229:
230: void ptoum(m,f,wf)
231: int m;
232: P f;
233: UM wf;
234: {
235: unsigned int r;
236: int i;
237: DCP dc;
238:
239: for ( i = UDEG(f); i >= 0; i-- )
240: COEF(wf)[i] = 0;
241:
242: for ( dc = DC(f); dc; dc = NEXT(dc) ) {
243: r = rem(NM((Q)COEF(dc)),m);
244: if ( r && (SGN((Q)COEF(dc)) < 0) )
245: r = m-r;
246: COEF(wf)[QTOS(DEG(dc))] = r;
247: }
248: degum(wf,UDEG(f));
249: }
250:
251: void umtop(v,w,f)
252: V v;
253: UM w;
254: P *f;
255: {
256: int *c;
257: DCP dc,dc0;
258: int i;
259: Q q;
260:
261: if ( DEG(w) < 0 )
262: *f = 0;
263: else if ( DEG(w) == 0 )
264: STOQ(COEF(w)[0],q), *f = (P)q;
265: else {
266: for ( i = DEG(w), c = COEF(w), dc0 = 0; i >= 0; i-- )
267: if ( c[i] ) {
268: NEXTDC(dc0,dc);
269: STOQ(i,DEG(dc));
270: STOQ(c[i],q), COEF(dc) = (P)q;
271: }
272: NEXT(dc) = 0;
273: MKP(v,dc0,*f);
274: }
275: }
276:
277: void ptoup(n,nr)
278: P n;
279: UP *nr;
280: {
281: DCP dc;
282: UP r;
283: int d;
284:
285: if ( !n )
286: *nr = 0;
287: else if ( OID(n) == O_N ) {
288: *nr = r = UPALLOC(0);
289: DEG(r) = 0; COEF(r)[0] = (Num)n;
290: } else {
291: d = UDEG(n);
292: up_var = VR(n);
293: *nr = r = UPALLOC(d); DEG(r) = d;
294: for ( dc = DC(n); dc; dc = NEXT(dc) ) {
295: COEF(r)[QTOS(DEG(dc))] = (Num)COEF(dc);
296: }
297: }
298: }
299:
300: void uptop(n,nr)
301: UP n;
302: P *nr;
303: {
304: int i;
305: DCP dc0,dc;
306:
307: if ( !n )
308: *nr = 0;
309: else if ( !DEG(n) )
310: *nr = (P)COEF(n)[0];
311: else {
312: for ( i = DEG(n), dc0 = 0; i >= 0; i-- )
313: if ( COEF(n)[i] ) {
314: NEXTDC(dc0,dc); STOQ(i,DEG(dc)); COEF(dc) = (P)COEF(n)[i];
315: }
316: if ( !up_var )
317: up_var = CO->v;
318: MKP(up_var,dc0,*nr);
319: }
320: }
321:
322: void ulmptoum(m,f,wf)
323: int m;
324: UP f;
325: UM wf;
326: {
327: int i,d;
328: LM *c;
329:
330: if ( !f )
331: wf->d = -1;
332: else {
333: wf->d = d = f->d;
334: c = (LM *)f->c;
335: for ( i = 0, d = f->d; i <= d; i++ )
336: COEF(wf)[i] = rem(c[i]->body,m);
337: }
338: }
339:
340: void objtobobj(base,p,rp)
341: int base;
342: Obj p;
343: Obj *rp;
344: {
345: if ( !p )
346: *rp = 0;
347: else
348: switch ( OID(p) ) {
349: case O_N:
350: numtobnum(base,(Num)p,(Num *)rp); break;
351: case O_P:
352: ptobp(base,(P)p,(P *)rp); break;
353: case O_LIST:
354: listtoblist(base,(LIST)p,(LIST *)rp); break;
355: case O_VECT:
356: vecttobvect(base,(VECT)p,(VECT *)rp); break;
357: case O_MAT:
358: mattobmat(base,(MAT)p,(MAT *)rp); break;
359: case O_STR:
360: *rp = p; break;
361: case O_COMP: default:
362: error("objtobobj : not implemented"); break;
363: }
364: }
365:
366: void bobjtoobj(base,p,rp)
367: int base;
368: Obj p;
369: Obj *rp;
370: {
371: if ( !p )
372: *rp = 0;
373: else
374: switch ( OID(p) ) {
375: case O_N:
376: bnumtonum(base,(Num)p,(Num *)rp); break;
377: case O_P:
378: bptop(base,(P)p,(P *)rp); break;
379: case O_LIST:
380: blisttolist(base,(LIST)p,(LIST *)rp); break;
381: case O_VECT:
382: bvecttovect(base,(VECT)p,(VECT *)rp); break;
383: case O_MAT:
384: bmattomat(base,(MAT)p,(MAT *)rp); break;
385: case O_STR:
386: *rp = p; break;
387: case O_COMP: default:
388: error("bobjtoobj : not implemented"); break;
389: }
390: }
391:
392: void numtobnum(base,p,rp)
393: int base;
394: Num p;
395: Num *rp;
396: {
397: N nm,dn,body;
398: Q q;
399: LM l;
400:
401: if ( !p )
402: *rp = 0;
403: else
404: switch ( NID(p) ) {
405: case N_Q:
406: ntobn(base,NM((Q)p),&nm);
407: if ( DN((Q)p) ) {
408: ntobn(base,DN((Q)p),&dn);
409: NDTOQ(nm,dn,SGN((Q)p),q);
410: } else
411: NTOQ(nm,SGN((Q)p),q);
412: *rp = (Num)q;
413: break;
414: case N_R:
415: *rp = p; break;
416: case N_LM:
417: ntobn(base,((LM)p)->body,&body);
418: MKLM(body,l); *rp = (Num)l;
419: break;
420: default:
421: error("numtobnum : not implemented"); break;
422: }
423: }
424:
425: void bnumtonum(base,p,rp)
426: int base;
427: Num p;
428: Num *rp;
429: {
430: N nm,dn,body;
431: Q q;
432: LM l;
433:
434: if ( !p )
435: *rp = 0;
436: else
437: switch ( NID(p) ) {
438: case N_Q:
439: bnton(base,NM((Q)p),&nm);
440: if ( DN((Q)p) ) {
441: bnton(base,DN((Q)p),&dn);
442: NDTOQ(nm,dn,SGN((Q)p),q);
443: } else
444: NTOQ(nm,SGN((Q)p),q);
445: *rp = (Num)q;
446: break;
447: case N_R:
448: *rp = p; break;
449: case N_LM:
450: bnton(base,((LM)p)->body,&body);
451: MKLM(body,l); *rp = (Num)l;
452: break;
453: default:
454: error("bnumtonum : not implemented"); break;
455: }
456: }
457:
458: void ptobp(base,p,rp)
459: int base;
460: P p;
461: P *rp;
462: {
463: DCP dcr0,dcr,dc;
464:
465: if ( !p )
466: *rp = p;
467: else {
468: for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
469: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc);
470: objtobobj(base,(Obj)COEF(dc),(Obj *)&COEF(dcr));
471: }
472: NEXT(dcr) = 0;
473: MKP(VR(p),dcr0,*rp);
474: }
475: }
476:
477: void bptop(base,p,rp)
478: int base;
479: P p;
480: P *rp;
481: {
482: DCP dcr0,dcr,dc;
483:
484: if ( !p )
485: *rp = p;
486: else {
487: for ( dcr0 = 0, dc = DC(p); dc; dc = NEXT(dc) ) {
488: NEXTDC(dcr0,dcr); DEG(dcr) = DEG(dc);
489: bobjtoobj(base,(Obj)COEF(dc),(Obj *)&COEF(dcr));
490: }
491: NEXT(dcr) = 0;
492: MKP(VR(p),dcr0,*rp);
493: }
494: }
495:
496: void listtoblist(base,p,rp)
497: int base;
498: LIST p;
499: LIST *rp;
500: {
501: NODE nr0,nr,n;
502:
503: if ( !p )
504: *rp = p;
505: else {
506: for ( nr0 = 0, n = BDY(p); n; n = NEXT(n) ) {
507: NEXTNODE(nr0,nr);
508: objtobobj(base,BDY(n),(Obj *)&BDY(nr));
509: }
510: NEXT(nr) = 0;
511: MKLIST(*rp,nr0);
512: }
513: }
514:
515: void blisttolist(base,p,rp)
516: int base;
517: LIST p;
518: LIST *rp;
519: {
520: NODE nr0,nr,n;
521:
522: if ( !p )
523: *rp = p;
524: else {
525: for ( nr0 = 0, n = BDY(p); n; n = NEXT(n) ) {
526: NEXTNODE(nr0,nr);
527: bobjtoobj(base,BDY(n),(Obj *)&BDY(nr));
528: }
529: NEXT(nr) = 0;
530: MKLIST(*rp,nr0);
531: }
532: }
533:
534: void vecttobvect(base,p,rp)
535: int base;
536: VECT p;
537: VECT *rp;
538: {
539: int i,l;
540: VECT r;
541:
542: if ( !p )
543: *rp = p;
544: else {
545: l = p->len;
546: MKVECT(r,l); *rp = r;
547: for ( i = 0; i < l; i++ )
548: objtobobj(base,p->body[i],(Obj *)&r->body[i]);
549: }
550: }
551:
552: void bvecttovect(base,p,rp)
553: int base;
554: VECT p;
555: VECT *rp;
556: {
557: int i,l;
558: VECT r;
559:
560: if ( !p )
561: *rp = p;
562: else {
563: l = p->len;
564: MKVECT(r,l); *rp = r;
565: for ( i = 0; i < l; i++ )
566: bobjtoobj(base,p->body[i],(Obj *)&r->body[i]);
567: }
568: }
569:
570: void mattobmat(base,p,rp)
571: int base;
572: MAT p;
573: MAT *rp;
574: {
575: int row,col,i,j;
576: MAT r;
577:
578: if ( !p )
579: *rp = p;
580: else {
581: row = p->row; col = p->col;
582: MKMAT(r,row,col); *rp = r;
583: for ( i = 0; i < row; i++ )
584: for ( j = 0; i < col; j++ )
585: objtobobj(base,p->body[i][j],(Obj *)&r->body[i][j]);
586: }
587: }
588:
589: void bmattomat(base,p,rp)
590: int base;
591: MAT p;
592: MAT *rp;
593: {
594: int row,col,i,j;
595: MAT r;
596:
597: if ( !p )
598: *rp = p;
599: else {
600: row = p->row; col = p->col;
601: MKMAT(r,row,col); *rp = r;
602: for ( i = 0; i < row; i++ )
603: for ( j = 0; i < col; j++ )
604: bobjtoobj(base,p->body[i][j],(Obj *)&r->body[i][j]);
605: }
606: }
607:
608: void n32ton27(g,rp)
609: N g;
610: N *rp;
611: {
612: int i,j,k,l,r,bits,words;
613: unsigned int t;
614: unsigned int *a,*b;
615: N z;
616:
617: l = PL(g); a = BD(g);
618: for ( i = 31, t = a[l-1]; !(t&(1<<i)); i-- );
619: bits = (l-1)*32+i+1; words = (bits+26)/27;
620: *rp = z = NALLOC(words); PL(z) = words;
621: bzero((char *)BD(z),words*sizeof(unsigned int));
622: for ( j = 0, b = BD(z); j < words; j++ ) {
623: k = (27*j)/32; r = (27*j)%32;
624: if ( r > 5 )
625: b[j] = (a[k]>>r)|(k==(l-1)?0:((a[k+1]&((1<<(r-5))-1))<<(32-r)));
626: else
627: b[j] = (a[k]>>r)&((1<<27)-1);
628: }
629: if ( !(r = bits%27) )
630: r = 27;
631: b[words-1] &= ((1<<r)-1);
632: }
633:
634: void n27ton32(a,rp)
635: N a;
636: N *rp;
637: {
638: int i,j,k,l,r,bits,words;
639: unsigned int t;
640: unsigned int *b,*c;
641: N z;
642:
643: l = PL(a); b = BD(a);
644: for ( i = 26, t = b[l-1]; !(t&(1<<i)); i-- );
645: bits = (l-1)*27+i+1; words = (bits+31)/32;
646: *rp = z = NALLOC(words); PL(z) = words;
647: bzero((char *)BD(z),words*sizeof(unsigned int));
648: for ( j = 0, c = BD(z); j < l; j++ ) {
649: k = (27*j)/32; r = (27*j)%32;
650: if ( r > 5 ) {
651: c[k] |= (b[j]&((1<<(32-r))-1))<<r;
652: if ( k+1 < words )
653: c[k+1] = (b[j]>>(32-r));
654: } else
655: c[k] |= (b[j]<<r);
656: }
657: }
658:
659: void mptoum(p,pr)
660: P p;
661: UM pr;
662: {
663: DCP dc;
664:
665: if ( !p )
666: DEG(pr) = -1;
667: else if ( NUM(p) ) {
668: DEG(pr) = 0; COEF(pr)[0] = CONT((MQ)p);
669: } else {
670: bzero((char *)pr,(int)((UDEG(p)+2)*sizeof(int)));
671: for ( dc = DC(p); dc; dc = NEXT(dc) )
672: COEF(pr)[QTOS(DEG(dc))] = CONT((MQ)COEF(dc));
673: degum(pr,UDEG(p));
674: }
675: }
676:
677: void umtomp(v,p,pr)
678: V v;
679: UM p;
680: P *pr;
681: {
682: DCP dc,dc0;
683: int i;
684: MQ q;
685:
686: if ( !p || (DEG(p) < 0) )
687: *pr = 0;
688: else if ( !DEG(p) )
689: STOMQ(COEF(p)[0],q), *pr = (P)q;
690: else {
691: for ( dc0 = 0, i = DEG(p); i >= 0; i-- )
692: if ( COEF(p)[i] ) {
693: NEXTDC(dc0,dc); STOQ(i,DEG(dc));
694: STOMQ(COEF(p)[i],q), COEF(dc) = (P)q;
695: }
696: NEXT(dc) = 0; MKP(v,dc0,*pr);
697: }
698: }