![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to N.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 N.c CVS log | Up to [local] / OpenXM_contrib2 / asir2000 / engine |
1.1 noro 1: /* $OpenXM: OpenXM/src/asir99/engine/N.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
2: #include "ca.h"
3: #include "base.h"
4:
5: void addn();
6: int subn();
7: void _bshiftn();
8: void dupn();
9:
10: #if defined(VISUAL) || defined(i386)
11: void addn(n1,n2,nr)
12: N n1,n2,*nr;
13: {
14: unsigned int *m1,*m2,*mr;
15: unsigned int c;
16: N r;
17: int i,d1,d2;
18: unsigned int tmp;
19:
20: if ( !n1 )
21: COPY(n2,*nr);
22: else if ( !n2 )
23: COPY(n1,*nr);
24: else {
25: if ( PL(n1) > PL(n2) ) {
26: d1 = PL(n1); d2 = PL(n2); m1 = BD(n1); m2 = BD(n2);
27: } else {
28: d1 = PL(n2); d2 = PL(n1); m1 = BD(n2); m2 = BD(n1);
29: }
30: *nr = r = NALLOC(d1 + 1); INITRC(r); mr = BD(r);
31:
32: #if defined(VISUAL)
33: __asm {
34: push esi
35: push edi
36: mov esi,m1
37: mov edi,m2
38: mov ebx,mr
39: mov ecx,d2
40: xor eax,eax
41: Lstart_addn:
42: mov eax,DWORD PTR [esi]
43: mov edx,DWORD PTR [edi]
44: adc eax,edx
45: mov DWORD PTR [ebx],eax
46: lea esi,DWORD PTR [esi+4]
47: lea edi,DWORD PTR [edi+4]
48: lea ebx,DWORD PTR [ebx+4]
49: dec ecx
50: jnz Lstart_addn
51: pop edi
52: pop esi
53: mov eax,0
54: adc eax,eax
55: mov c,eax
56: }
57: #else
58: asm volatile("\
59: movl %1,%%esi;\
60: movl %2,%%edi;\
61: movl %3,%%ebx;\
62: movl %4,%%ecx;\
63: testl %%eax,%%eax;\
64: Lstart_addn:;\
65: movl (%%esi),%%eax;\
66: movl (%%edi),%%edx;\
67: adcl %%edx,%%eax;\
68: movl %%eax,(%%ebx);\
69: leal 4(%%esi),%%esi;\
70: leal 4(%%edi),%%edi;\
71: leal 4(%%ebx),%%ebx;\
72: decl %%ecx;\
73: jnz Lstart_addn;\
74: movl $0,%%eax;\
75: adcl %%eax,%%eax;\
76: movl %%eax,%0"\
77: :"=m"(c)\
78: :"m"(m1),"m"(m2),"m"(mr),"m"(d2)\
79: :"eax","ebx","ecx","edx","esi","edi");
80: #endif
81: for ( i = d2, m1 += d2, mr += d2; (i < d1) && c ; i++ ) {
82: tmp = *m1++ + c;
83: c = tmp < c ? 1 : 0;
84: *mr++ = tmp;
85: }
86: for ( ; i < d1; i++ )
87: *mr++ = *m1++;
88: *mr = c;
89: PL(r) = (c?d1+1:d1);
90: }
91: }
92:
93: int subn(n1,n2,nr)
94: N n1,n2,*nr;
95: {
96: N r;
97: unsigned int *m1,*m2,*mr,br;
98: unsigned int tmp,t;
99: int d1,d2,sgn,i;
100:
101: if ( !n1 ) {
102: if ( n2 ) {
103: COPY(n2,*nr);
104: return -1;
105: } else {
106: *nr = 0;
107: return 0;
108: }
109: } else if ( !n2 ) {
110: COPY(n1,*nr);
111: return 1;
112: } else {
113: d1 = PL(n1); d2 = PL(n2);
114: m1 = BD(n1); m2 = BD(n2);
115: if ( (d1 = PL(n1)) > (d2 = PL(n2)) )
116: sgn = 1;
117: else if ( d1 < d2 ) {
118: d1 = PL(n2); d2 = PL(n1);
119: m1 = BD(n2); m2 = BD(n1);
120: sgn = -1;
121: } else {
122: for ( i = d1-1; i >= 0 && m1[i] == m2[i]; i-- );
123: if ( i < 0 ) {
124: *nr = 0;
125: return 0;
126: }
127: d1 = d2 = i+1;
128: if ( m1[i] > m2[i] )
129: sgn = 1;
130: else {
131: m1 = BD(n2); m2 = BD(n1);
132: sgn = -1;
133: }
134: }
135: *nr = r = NALLOC(d1); INITRC(r); mr = BD(r);
136:
137: #if defined(VISUAL)
138: __asm {
139: push esi
140: push edi
141: mov esi,m1
142: mov edi,m2
143: mov ebx,mr
144: mov ecx,d2
145: xor eax,eax
146: Lstart_subn:
147: mov eax,DWORD PTR [esi]
148: mov edx,DWORD PTR [edi]
149: sbb eax,edx
150: mov DWORD PTR [ebx],eax
151: lea esi,DWORD PTR [esi+4]
152: lea edi,DWORD PTR [edi+4]
153: lea ebx,DWORD PTR [ebx+4]
154: dec ecx
155: jnz Lstart_subn
156: pop edi
157: pop esi
158: mov eax,0
159: adc eax,eax
160: mov br,eax
161: }
162: #else
163: asm volatile("\
164: movl %1,%%esi;\
165: movl %2,%%edi;\
166: movl %3,%%ebx;\
167: movl %4,%%ecx;\
168: testl %%eax,%%eax;\
169: Lstart_subn:;\
170: movl (%%esi),%%eax;\
171: movl (%%edi),%%edx;\
172: sbbl %%edx,%%eax;\
173: movl %%eax,(%%ebx);\
174: leal 4(%%esi),%%esi;\
175: leal 4(%%edi),%%edi;\
176: leal 4(%%ebx),%%ebx;\
177: decl %%ecx;\
178: jnz Lstart_subn;\
179: movl $0,%%eax;\
180: adcl %%eax,%%eax;\
181: movl %%eax,%0"\
182: :"=m"(br)\
183: :"m"(m1),"m"(m2),"m"(mr),"m"(d2)\
184: :"eax","ebx","ecx","edx","esi","edi");
185: #endif
186: for ( i = d2, m1 += d2, mr += d2; (i < d1) && br; i++ ) {
187: t = *m1++;
188: tmp = t - br;
189: br = tmp > t ? 1 : 0;
190: *mr++ = tmp;
191: }
192: for ( ; i < d1; i++ )
193: *mr++ = *m1++;
194: for ( i = d1-1, mr--; i >= 0 && !*mr--; i-- );
195: PL(r) = i + 1;
196: return sgn;
197: }
198: }
199:
200: void _addn(n1,n2,nr)
201: N n1,n2,nr;
202: {
203: unsigned int *m1,*m2,*mr;
204: unsigned int c;
205: int i,d1,d2;
206: unsigned int tmp;
207:
208: if ( !n1 || !PL(n1) )
209: dupn(n2,nr);
210: else if ( !n2 || !PL(n2) )
211: dupn(n1,nr);
212: else {
213: if ( PL(n1) > PL(n2) ) {
214: d1 = PL(n1); d2 = PL(n2); m1 = BD(n1); m2 = BD(n2);
215: } else {
216: d1 = PL(n2); d2 = PL(n1); m1 = BD(n2); m2 = BD(n1);
217: }
218: mr = BD(nr);
219:
220: #if defined(VISUAL)
221: __asm {
222: push esi
223: push edi
224: mov esi,m1
225: mov edi,m2
226: mov ebx,mr
227: mov ecx,d2
228: xor eax,eax
229: Lstart__addn:
230: mov eax,DWORD PTR [esi]
231: mov edx,DWORD PTR [edi]
232: adc eax,edx
233: mov DWORD PTR [ebx],eax
234: lea esi,DWORD PTR [esi+4]
235: lea edi,DWORD PTR [edi+4]
236: lea ebx,DWORD PTR [ebx+4]
237: dec ecx
238: jnz Lstart__addn
239: pop edi
240: pop esi
241: mov eax,0
242: adc eax,eax
243: mov c,eax
244: }
245: #else
246: asm volatile("\
247: movl %1,%%esi;\
248: movl %2,%%edi;\
249: movl %3,%%ebx;\
250: movl %4,%%ecx;\
251: testl %%eax,%%eax;\
252: Lstart__addn:;\
253: movl (%%esi),%%eax;\
254: movl (%%edi),%%edx;\
255: adcl %%edx,%%eax;\
256: movl %%eax,(%%ebx);\
257: leal 4(%%esi),%%esi;\
258: leal 4(%%edi),%%edi;\
259: leal 4(%%ebx),%%ebx;\
260: decl %%ecx;\
261: jnz Lstart__addn;\
262: movl $0,%%eax;\
263: adcl %%eax,%%eax;\
264: movl %%eax,%0"\
265: :"=m"(c)\
266: :"m"(m1),"m"(m2),"m"(mr),"m"(d2)\
267: :"eax","ebx","ecx","edx","esi","edi");
268: #endif
269: for ( i = d2, m1 += d2, mr += d2; (i < d1) && c ; i++ ) {
270: tmp = *m1++ + c;
271: c = tmp < c ? 1 : 0;
272: *mr++ = tmp;
273: }
274: for ( ; i < d1; i++ )
275: *mr++ = *m1++;
276: *mr = c;
277: PL(nr) = (c?d1+1:d1);
278: }
279: }
280:
281: int _subn(n1,n2,nr)
282: N n1,n2,nr;
283: {
284: unsigned int *m1,*m2,*mr,br;
285: unsigned int tmp,t;
286: int d1,d2,sgn,i;
287:
288: if ( !n1 || !PL(n1) ) {
289: if ( n2 && PL(n2) ) {
290: dupn(n2,nr);
291: return -1;
292: } else {
293: PL(nr) = 0;
294: return 0;
295: }
296: } else if ( !n2 || !PL(n2) ) {
297: dupn(n1,nr);
298: return 1;
299: } else {
300: d1 = PL(n1); d2 = PL(n2);
301: m1 = BD(n1); m2 = BD(n2);
302: if ( (d1 = PL(n1)) > (d2 = PL(n2)) )
303: sgn = 1;
304: else if ( d1 < d2 ) {
305: d1 = PL(n2); d2 = PL(n1);
306: m1 = BD(n2); m2 = BD(n1);
307: sgn = -1;
308: } else {
309: for ( i = d1-1; i >= 0 && m1[i] == m2[i]; i-- );
310: if ( i < 0 ) {
311: PL(nr) = 0;
312: return 0;
313: }
314: d1 = d2 = i+1;
315: if ( m1[i] > m2[i] )
316: sgn = 1;
317: else {
318: m1 = BD(n2); m2 = BD(n1);
319: sgn = -1;
320: }
321: }
322: mr = BD(nr);
323:
324: #if defined(VISUAL)
325: __asm {
326: push esi
327: push edi
328: mov esi,m1
329: mov edi,m2
330: mov ebx,mr
331: mov ecx,d2
332: xor eax,eax
333: Lstart__subn:
334: mov eax,DWORD PTR [esi]
335: mov edx,DWORD PTR [edi]
336: sbb eax,edx
337: mov DWORD PTR [ebx],eax
338: lea esi,DWORD PTR [esi+4]
339: lea edi,DWORD PTR [edi+4]
340: lea ebx,DWORD PTR [ebx+4]
341: dec ecx
342: jnz Lstart__subn
343: pop edi
344: pop esi
345: mov eax,0
346: adc eax,eax
347: mov br,eax
348: }
349: #else
350: asm volatile("\
351: movl %1,%%esi;\
352: movl %2,%%edi;\
353: movl %3,%%ebx;\
354: movl %4,%%ecx;\
355: testl %%eax,%%eax;\
356: Lstart__subn:;\
357: movl (%%esi),%%eax;\
358: movl (%%edi),%%edx;\
359: sbbl %%edx,%%eax;\
360: movl %%eax,(%%ebx);\
361: leal 4(%%esi),%%esi;\
362: leal 4(%%edi),%%edi;\
363: leal 4(%%ebx),%%ebx;\
364: decl %%ecx;\
365: jnz Lstart__subn;\
366: movl $0,%%eax;\
367: adcl %%eax,%%eax;\
368: movl %%eax,%0"\
369: :"=m"(br)\
370: :"m"(m1),"m"(m2),"m"(mr),"m"(d2)\
371: :"eax","ebx","ecx","edx","esi","edi");
372: #endif
373: for ( i = d2, m1 += d2, mr += d2; (i < d1) && br; i++ ) {
374: t = *m1++;
375: tmp = t - br;
376: br = tmp > t ? 1 : 0;
377: *mr++ = tmp;
378: }
379: for ( ; i < d1; i++ )
380: *mr++ = *m1++;
381: for ( i = d1-1, mr--; i >= 0 && !*mr--; i-- );
382: PL(nr) = i + 1;
383: return sgn;
384: }
385: }
386: #else
387:
388: void addn(n1,n2,nr)
389: N n1,n2,*nr;
390: {
391: unsigned int *m1,*m2,*mr,i,c;
392: N r;
393: int d1,d2;
394: unsigned int tmp;
395:
396: if ( !n1 )
397: COPY(n2,*nr);
398: else if ( !n2 )
399: COPY(n1,*nr);
400: else {
401: if ( PL(n1) > PL(n2) ) {
402: d1 = PL(n1); d2 = PL(n2); m1 = BD(n1); m2 = BD(n2);
403: } else {
404: d1 = PL(n2); d2 = PL(n1); m1 = BD(n2); m2 = BD(n1);
405: }
406: *nr = r = NALLOC(d1 + 1); INITRC(r);
407: for ( i = 0, c = 0, mr = BD(r); i < d2; i++, m1++, m2++, mr++ ) {
408: tmp = *m1 + *m2;
409: if ( tmp < *m1 ) {
410: tmp += c;
411: c = 1;
412: } else {
413: tmp += c;
414: c = tmp < c ? 1 : 0;
415: }
416: *mr = tmp;
417: }
418: for ( ; (i < d1) && c ; i++, m1++, mr++ ) {
419: tmp = *m1 + c;
420: c = tmp < c ? 1 : 0;
421: *mr = tmp;
422: }
423: for ( ; i < d1; i++ )
424: *mr++ = *m1++;
425: *mr = c;
426: PL(r) = (c?d1+1:d1);
427: }
428: }
429:
430: int subn(n1,n2,nr)
431: N n1,n2,*nr;
432: {
433: N r;
434: unsigned int *m1,*m2,*mr,i,br;
435: L tmp;
436: int d1,d2,nz,sgn;
437:
438: if ( !n1 ) {
439: if ( n2 ) {
440: COPY(n2,*nr);
441: return -1;
442: } else {
443: *nr = 0;
444: return 0;
445: }
446: } else if ( !n2 ) {
447: COPY(n1,*nr);
448: return 1;
449: } else {
450: switch ( cmpn(n1,n2) ) {
451: case 1:
452: d1 = PL(n1); d2 = PL(n2); m1 = BD(n1); m2 = BD(n2);
453: sgn = 1; break;
454: case -1:
455: d1 = PL(n2); d2 = PL(n1); m1 = BD(n2); m2 = BD(n1);
456: sgn = -1; break;
457: case 0:
458: default:
459: *nr = 0; return ( 0 ); break;
460: }
461: *nr = r = NALLOC(d1); INITRC(r);
462: for ( i = 0, br = 0, nz = -1, mr = BD(r);
463: i < d2; i++ ) {
464: if ( (tmp = (L)*m1++ - (L)*m2++ - (L)br) && ( tmp > -LBASE ) )
465: nz = i;
466: if ( tmp < 0 ) {
467: br = 1; *mr++ = (unsigned int)(tmp + LBASE);
468: } else {
469: br = 0; *mr++ = (unsigned int)tmp;
470: }
471: }
472: for ( ; (i < d1) && br; i++ ) {
473: if ( (tmp = (L)*m1++ - (L)br) && ( tmp > -LBASE ) )
474: nz = i;
475: if ( tmp < 0 ) {
476: br = 1; *mr++ = (unsigned int)(tmp + LBASE);
477: } else {
478: br = 0; *mr++ = (unsigned int)tmp;
479: }
480: }
481: for ( ; i < d1; i++ )
482: if ( *mr++ = *m1++ )
483: nz = i;
484: PL(r) = nz + 1;
485: return sgn;
486: }
487: }
488:
489: void _addn(n1,n2,nr)
490: N n1,n2,nr;
491: {
492: unsigned int *m1,*m2,*mr,i,c;
493: int d1,d2;
494: unsigned int tmp;
495:
496: if ( !n1 || !PL(n1) )
497: dupn(n2,nr);
498: else if ( !n2 || !PL(n2) )
499: dupn(n1,nr);
500: else {
501: if ( PL(n1) > PL(n2) ) {
502: d1 = PL(n1); d2 = PL(n2); m1 = BD(n1); m2 = BD(n2);
503: } else {
504: d1 = PL(n2); d2 = PL(n1); m1 = BD(n2); m2 = BD(n1);
505: }
506: for ( i = 0, c = 0, mr = BD(nr); i < d2; i++, m1++, m2++, mr++ ) {
507: tmp = *m1 + *m2;
508: if ( tmp < *m1 ) {
509: tmp += c;
510: c = 1;
511: } else {
512: tmp += c;
513: c = tmp < c ? 1 : 0;
514: }
515: *mr = tmp;
516: }
517: for ( ; (i < d1) && c ; i++, m1++, mr++ ) {
518: tmp = *m1 + c;
519: c = tmp < c ? 1 : 0;
520: *mr = tmp;
521: }
522: for ( ; i < d1; i++ )
523: *mr++ = *m1++;
524: *mr = c;
525: PL(nr) = (c?d1+1:d1);
526: }
527: }
528:
529: int _subn(n1,n2,nr)
530: N n1,n2,nr;
531: {
532: N r;
533: unsigned int *m1,*m2,*mr,i,br;
534: L tmp;
535: int d1,d2,nz,sgn;
536:
537: if ( !n1 || !PL(n1) ) {
538: if ( n2 && PL(n2) ) {
539: dupn(n2,nr);
540: return -1;
541: } else {
542: PL(nr) = 0;
543: return 0;
544: }
545: } else if ( !n2 || !PL(n2) ) {
546: dupn(n1,nr);
547: return 1;
548: } else {
549: switch ( cmpn(n1,n2) ) {
550: case 1:
551: d1 = PL(n1); d2 = PL(n2); m1 = BD(n1); m2 = BD(n2);
552: sgn = 1; break;
553: case -1:
554: d1 = PL(n2); d2 = PL(n1); m1 = BD(n2); m2 = BD(n1);
555: sgn = -1; break;
556: case 0:
557: default:
558: PL(nr) = 0; return ( 0 ); break;
559: }
560: for ( i = 0, br = 0, nz = -1, mr = BD(nr);
561: i < d2; i++ ) {
562: if ( (tmp = (L)*m1++ - (L)*m2++ - (L)br) && ( tmp > -LBASE ) )
563: nz = i;
564: if ( tmp < 0 ) {
565: br = 1; *mr++ = (unsigned int)(tmp + LBASE);
566: } else {
567: br = 0; *mr++ = (unsigned int)tmp;
568: }
569: }
570: for ( ; (i < d1) && br; i++ ) {
571: if ( (tmp = (L)*m1++ - (L)br) && ( tmp > -LBASE ) )
572: nz = i;
573: if ( tmp < 0 ) {
574: br = 1; *mr++ = (unsigned int)(tmp + LBASE);
575: } else {
576: br = 0; *mr++ = (unsigned int)tmp;
577: }
578: }
579: for ( ; i < d1; i++ )
580: if ( *mr++ = *m1++ )
581: nz = i;
582: PL(nr) = nz + 1;
583: return sgn;
584: }
585: }
586: #endif
587:
588: /* a2 += a1; n2 >= n1 */
589:
590: void addarray_to(a1,n1,a2,n2)
591: unsigned int *a1;
592: int n1;
593: unsigned int *a2;
594: int n2;
595: {
596: int i;
597: unsigned int c,tmp;
598:
599: for ( i = 0, c = 0; i < n2; i++, a1++, a2++ ) {
600: tmp = *a1 + *a2;
601: if ( tmp < *a1 ) {
602: tmp += c;
603: c = 1;
604: } else {
605: tmp += c;
606: c = tmp < c ? 1 : 0;
607: }
608: *a2 = tmp;
609: }
610: for ( ; (i < n2) && c ; i++, a2++ ) {
611: tmp = *a2 + c;
612: c = tmp < c ? 1 : 0;
613: *a2 = tmp;
614: }
615: if ( i == n2 && c )
616: *a2 = c;
617: }
618:
619: void pwrn(n,e,nr)
620: N n,*nr;
621: int e;
622: {
623: N nw,nw1;
624:
625: if ( e == 1 ) {
626: COPY(n,*nr);
627: return;
628: }
629: pwrn(n,e / 2,&nw);
630: if ( e % 2 == 0 )
631: kmuln(nw,nw,nr);
632: else {
633: kmuln(nw,nw,&nw1); kmuln(nw1,n,nr); FREEN(nw1);
634: }
635: FREEN(nw);
636: }
637:
638: extern int igcd_algorithm;
639:
640: void gcdEuclidn();
641:
642: void gcdn(n1,n2,nr)
643: N n1,n2,*nr;
644: {
645: N m1,m2,g;
646:
647: if ( !igcd_algorithm )
648: gcdEuclidn(n1,n2,nr);
649: else {
650: if ( !n1 )
651: *nr = n2;
652: else if ( !n2 )
653: *nr = n1;
654: else {
655: n32ton27(n1,&m1);
656: n32ton27(n2,&m2);
657: gcdBinary_27n(m1,m2,&g);
658: n27ton32(g,nr);
659: }
660: }
661: }
662:
663: void gcdEuclidn(n1,n2,nr)
664: N n1,n2,*nr;
665: {
666: N m1,m2,q,r;
667: unsigned int i1,i2,ir;
668:
669: if ( !n1 )
670: COPY(n2,*nr);
671: else if ( !n2 )
672: COPY(n1,*nr);
673: else {
674: if ( PL(n1) > PL(n2) ) {
675: COPY(n1,m1); COPY(n2,m2);
676: } else {
677: COPY(n1,m2); COPY(n2,m1);
678: }
679: while ( PL(m1) > 1 ) {
680: divn(m1,m2,&q,&r); FREEN(m1); FREEN(q);
681: if ( !r ) {
682: *nr = m2;
683: return;
684: } else {
685: m1 = m2; m2 = r;
686: }
687: }
688: for ( i1 = BD(m1)[0], i2 = BD(m2)[0]; ir = i1 % i2; i1 = i2, i2 = ir );
689: if ( i2 == 1 )
690: COPY(ONEN,*nr);
691: else {
692: *nr = r = NALLOC(1); INITRC(r); PL(r) = 1; BD(r)[0] = i2;
693: }
694: }
695: }
696:
697: int cmpn(n1,n2)
698: N n1,n2;
699: {
700: int i;
701: unsigned int *m1,*m2;
702:
703: if ( !n1 )
704: if ( !n2 )
705: return 0;
706: else
707: return -1;
708: else if ( !n2 )
709: return 1;
710: else if ( PL(n1) > PL(n2) )
711: return 1;
712: else if ( PL(n1) < PL(n2) )
713: return -1;
714: else {
715: for ( i = PL(n1)-1, m1 = BD(n1)+i, m2 = BD(n2)+i;
716: i >= 0; i--, m1--, m2-- )
717: if ( *m1 > *m2 )
718: return 1;
719: else if ( *m1 < *m2 )
720: return -1;
721: return 0;
722: }
723: }
724:
725: void bshiftn(n,b,r)
726: N n;
727: int b;
728: N *r;
729: {
730: int w,l,nl,i,j;
731: N z;
732: unsigned int msw;
733: unsigned int *p,*pz;
734:
735: if ( b == 0 ) {
736: COPY(n,*r); return;
737: }
738: if ( b > 0 ) { /* >> */
739: w = b / BSH; l = PL(n)-w;
740: if ( l <= 0 ) {
741: *r = 0; return;
742: }
743: b %= BSH; p = BD(n)+w;
744: if ( !b ) {
745: *r = z = NALLOC(l); INITRC(z); PL(z) = l;
746: bcopy(p,BD(z),l*sizeof(unsigned int));
747: return;
748: }
749: msw = p[l-1];
750: for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
751: if ( b >= i ) {
752: l--;
753: if ( !l ) {
754: *r = 0; return;
755: }
756: *r = z = NALLOC(l); INITRC(z); PL(z) = l; pz = BD(z);
757: for ( j = 0; j < l; j++, p++ )
758: *pz++ = (*(p+1)<<(BSH-b))|(*p>>b);
759: } else {
760: *r = z = NALLOC(l); INITRC(z); PL(z) = l; pz = BD(z);
761: for ( j = 1; j < l; j++, p++ )
762: *pz++ = (*(p+1)<<(BSH-b))|(*p>>b);
763: *pz = *p>>b;
764: }
765: } else { /* << */
766: b = -b;
767: w = b / BSH; b %= BSH; l = PL(n); p = BD(n);
768: if ( !b ) {
769: nl = l+w; *r = z = NALLOC(nl); INITRC(z); PL(z) = nl;
770: bzero((char *)BD(z),w*sizeof(unsigned int));
771: bcopy(p,BD(z)+w,l*sizeof(unsigned int));
772: return;
773: }
774: msw = p[l-1];
775: for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
776: if ( b + i > BSH ) {
777: nl = l+w+1;
778: *r = z = NALLOC(nl); INITRC(z); PL(z) = nl; pz = BD(z)+w;
779: bzero((char *)BD(z),w*sizeof(unsigned int));
780: *pz++ = *p++<<b;
781: for ( j = 1; j < l; j++, p++ )
782: *pz++ = (*p<<b)|(*(p-1)>>(BSH-b));
783: *pz = *(p-1)>>(BSH-b);
784: } else {
785: nl = l+w;
786: *r = z = NALLOC(nl); INITRC(z); PL(z) = nl; pz = BD(z)+w;
787: bzero((char *)BD(z),w*sizeof(unsigned int));
788: *pz++ = *p++<<b;
789: for ( j = 1; j < l; j++, p++ )
790: *pz++ = (*p<<b)|(*(p-1)>>(BSH-b));
791: }
792: }
793: }
794:
795: #if 0
796: void _bshiftn(n,b,z)
797: N n;
798: int b;
799: N z;
800: {
801: int w,l,nl,i,j;
802: unsigned int msw;
803: unsigned int *p,*pz;
804:
805: if ( b == 0 ) {
806: copyn(n,PL(n),BD(z)); PL(z) = PL(n); return;
807: }
808: if ( b > 0 ) { /* >> */
809: w = b / BSH; l = PL(n)-w;
810: if ( l <= 0 ) {
811: PL(z) = 0; return;
812: }
813: b %= BSH; p = BD(n)+w;
814: if ( !b ) {
815: PL(z) = l;
816: bcopy(p,BD(z),l*sizeof(unsigned int));
817: return;
818: }
819: msw = p[l-1];
820: for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
821: if ( b >= i ) {
822: l--;
823: if ( !l ) {
824: PL(z) = 0; return;
825: }
826: PL(z) = l; pz = BD(z);
827: for ( j = 0; j < l; j++, p++ )
828: *pz++ = (*(p+1)<<(BSH-b))|(*p>>b);
829: } else {
830: PL(z) = l; pz = BD(z);
831: for ( j = 1; j < l; j++, p++ )
832: *pz++ = (*(p+1)<<(BSH-b))|(*p>>b);
833: *pz = *p>>b;
834: }
835: } else { /* << */
836: b = -b;
837: w = b / BSH; b %= BSH; l = PL(n); p = BD(n);
838: if ( !b ) {
839: nl = l+w; PL(z) = nl;
840: bzero((char *)BD(z),w*sizeof(unsigned int));
841: bcopy(p,BD(z)+w,l*sizeof(unsigned int));
842: return;
843: }
844: msw = p[l-1];
845: for ( i = BSH-1; !(msw&(((unsigned int)1)<<i)); i-- ); i++;
846: if ( b + i > BSH ) {
847: nl = l+w+1;
848: PL(z) = nl; pz = BD(z)+w;
849: bzero((char *)BD(z),w*sizeof(unsigned int));
850: *pz++ = *p++<<b;
851: for ( j = 1; j < l; j++, p++ )
852: *pz++ = (*p<<b)|(*(p-1)>>(BSH-b));
853: *pz = *(p-1)>>(BSH-b);
854: } else {
855: nl = l+w;
856: PL(z) = nl; pz = BD(z)+w;
857: bzero((char *)BD(z),w*sizeof(unsigned int));
858: *pz++ = *p++<<b;
859: for ( j = 1; j < l; j++, p++ )
860: *pz++ = (*p<<b)|(*(p-1)>>(BSH-b));
861: }
862: }
863: }
864: #endif
865:
866: void shiftn(n,w,r)
867: N n;
868: int w;
869: N *r;
870: {
871: int l,nl;
872: N z;
873:
874: if ( w == 0 )
875: COPY(n,*r);
876: else if ( w > 0 ) { /* >> */
877: l = PL(n)-w;
878: if ( l <= 0 )
879: *r = 0;
880: else {
881: *r = z = NALLOC(l); INITRC(z); PL(z) = l;
882: bcopy(BD(n)+w,BD(z),l*sizeof(unsigned int));
883: }
884: } else { /* << */
885: w = -w;
886: l = PL(n); nl = l+w;
887: *r = z = NALLOC(nl); INITRC(z); PL(z) = nl;
888: bzero((char *)BD(z),w*sizeof(unsigned int));
889: bcopy(BD(n),BD(z)+w,l*sizeof(unsigned int));
890: }
891: }
892:
893: void randomn(bits,r)
894: int bits;
895: N *r;
896: {
897: int l,i;
898: unsigned int *tb;
899: N t;
900:
901: l = (bits+31)>>5; /* word length */
902: *r = t = NALLOC(l);
903: tb = BD(t);
904: for ( i = 0; i < l; i++ )
905: tb[i] = mt_genrand();
906: if ( bits&31 )
907: tb[l-1] &= (1<<(bits&31))-1;
908: for ( i = l-1; i >= 0 && !tb[i]; i-- );
909: if ( i < 0 )
910: *r = 0;
911: else
912: PL(t) = i+1;
913: }
914:
915: void freen(n)
916: N n;
917: {
918: if ( n && (n != ONEN) )
919: free(n);
920: }
921:
922: int n_bits(n)
923: N n;
924: {
925: unsigned int l,i,t;
926:
927: if ( !n )
928: return 0;
929: l = PL(n); t = BD(n)[l-1];
930: for ( i = 0; t; t>>=1, i++);
931: return i + (l-1)*BSH;
932: }