Annotation of OpenXM_contrib2/asir2000/builtin/pdiv.c, Revision 1.7
1.3 noro 1: /*
2: * Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED
3: * All rights reserved.
4: *
5: * FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited,
6: * non-exclusive and royalty-free license to use, copy, modify and
7: * redistribute, solely for non-commercial and non-profit purposes, the
8: * computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and
9: * conditions of this Agreement. For the avoidance of doubt, you acquire
10: * only a limited right to use the SOFTWARE hereunder, and FLL or any
11: * third party developer retains all rights, including but not limited to
12: * copyrights, in and to the SOFTWARE.
13: *
14: * (1) FLL does not grant you a license in any way for commercial
15: * purposes. You may use the SOFTWARE only for non-commercial and
16: * non-profit purposes only, such as academic, research and internal
17: * business use.
18: * (2) The SOFTWARE is protected by the Copyright Law of Japan and
19: * international copyright treaties. If you make copies of the SOFTWARE,
20: * with or without modification, as permitted hereunder, you shall affix
21: * to all such copies of the SOFTWARE the above copyright notice.
22: * (3) An explicit reference to this SOFTWARE and its copyright owner
23: * shall be made on your publication or presentation in any form of the
24: * results obtained by use of the SOFTWARE.
25: * (4) In the event that you modify the SOFTWARE, you shall notify FLL by
1.4 noro 26: * e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification
1.3 noro 27: * for such modification or the source code of the modified part of the
28: * SOFTWARE.
29: *
30: * THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL
31: * MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND
32: * EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS
33: * FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES'
34: * RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY
35: * MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY.
36: * UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT,
37: * OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY
38: * DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL
39: * DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES
40: * ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES
41: * FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY
42: * DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF
43: * SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART
44: * OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY
45: * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
46: * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
47: *
1.7 ! noro 48: * $OpenXM: OpenXM_contrib2/asir2000/builtin/pdiv.c,v 1.6 2001/03/29 09:49:56 noro Exp $
1.3 noro 49: */
1.1 noro 50: #include "ca.h"
51: #include "parse.h"
52:
53: void Psdiv(), Psrem(), Ptdiv(), Psqr(), Pinva_mod();
1.7 ! noro 54: void Psdiv_gf2n(), Psrem_gf2n(), Pgcd_gf2n();
1.1 noro 55: void Psdivm(), Psremm(), Psqrm();
56: void Psrem_mod();
57: void Pugcd();
58: void Purem();
59: void Pudiv();
60:
61: struct ftab pdiv_tab[] = {
62: {"sdiv",Psdiv,-3},
63: {"srem",Psrem,-3},
64: {"sdiv_gf2n",Psdiv_gf2n,2},
65: {"srem_gf2n",Psrem_gf2n,2},
1.7 ! noro 66: {"gcd_gf2n",Pgcd_gf2n,2},
1.1 noro 67: {"sqr",Psqr,-3},
1.5 noro 68: {"tdiv",Ptdiv,-3},
1.1 noro 69: {"udiv",Pudiv,2},
70: {"sdivm",Psdivm,-4},
71: {"sremm",Psremm,-4},
72: {"sqrm",Psqrm,-4},
73: {"inva_mod",Pinva_mod,3},
74: {"srem_mod",Psrem_mod,3},
75: {"ugcd",Pugcd,2},
76: {"urem",Purem,2},
77: {0,0,0},
78: };
79:
80: void Psdiv(arg,rp)
81: NODE arg;
82: Obj *rp;
83: {
84: P q,r,dnd,dnd1,dvr,dvr1;
85: V v;
86: VL vl;
87:
88: asir_assert(ARG0(arg),O_P,"sdiv");
89: asir_assert(ARG1(arg),O_P,"sdiv");
90: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
91: if ( argc(arg) == 3 ) {
92: v = VR((P)ARG2(arg));
93: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
94: reordvar(CO,v,&vl);
95: divsrp(vl,dnd1,dvr1,&q,&r);
96: restore_mvar(CO,q,v,(P *)rp);
97: } else
98: divsrp(CO,dnd,dvr,(P *)rp,&r);
99: }
100:
101: void Psrem(arg,rp)
102: NODE arg;
103: Obj *rp;
104: {
105: P q,r,dnd,dnd1,dvr,dvr1;
106: V v;
107: VL vl;
108:
109: asir_assert(ARG0(arg),O_P,"srem");
110: asir_assert(ARG1(arg),O_P,"srem");
111: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
112: if ( argc(arg) == 3 ) {
113: v = VR((P)ARG2(arg));
114: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
115: reordvar(CO,v,&vl);
116: divsrp(vl,dnd1,dvr1,&q,&r);
117: restore_mvar(CO,r,v,(P *)rp);
118: } else
119: divsrp(CO,dnd,dvr,&q,(P *)rp);
120: }
121:
122: void Psqr(arg,rp)
123: NODE arg;
124: LIST *rp;
125: {
126: P q,q1,r,r1,dnd,dnd1,dvr,dvr1;
127: NODE n,tn;
128: V v;
129: VL vl;
130:
131: asir_assert(ARG0(arg),O_P,"sqr");
132: asir_assert(ARG1(arg),O_P,"sqr");
133: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
134: if ( argc(arg) == 3 ) {
135: v = VR((P)ARG2(arg));
136: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
137: reordvar(CO,v,&vl);
138: divsrp(vl,dnd1,dvr1,&q1,&r1);
139: restore_mvar(CO,q1,v,&q); restore_mvar(CO,r1,v,&r);
140: } else
141: divsrp(CO,dnd,dvr,&q,&r);
142: MKNODE(tn,r,0); MKNODE(n,q,tn); MKLIST(*rp,n);
143: }
144:
145: void Psdiv_gf2n(arg,rp)
146: NODE arg;
147: GF2N *rp;
148: {
149: GF2N dnd,dvr;
150: UP2 q,r;
151:
152: dnd = (GF2N)ARG0(arg); dvr = (GF2N)ARG1(arg);
153: if ( !dvr )
154: error("sdiv_gf2n : division by 0");
155: else if ( !dnd )
156: *rp = 0;
157: else {
158: qrup2(dnd->body,dvr->body,&q,&r);
159: MKGF2N(q,*rp);
160: }
161: }
162:
163: void Psrem_gf2n(arg,rp)
164: NODE arg;
165: GF2N *rp;
166: {
167: GF2N dnd,dvr;
168: UP2 q,r;
169:
170: dnd = (GF2N)ARG0(arg); dvr = (GF2N)ARG1(arg);
171: if ( !dvr )
172: error("srem_gf2n : division by 0");
173: else if ( !dnd )
174: *rp = 0;
175: else {
176: qrup2(dnd->body,dvr->body,&q,&r);
177: MKGF2N(r,*rp);
1.7 ! noro 178: }
! 179: }
! 180:
! 181: void Pgcd_gf2n(arg,rp)
! 182: NODE arg;
! 183: GF2N *rp;
! 184: {
! 185: GF2N p1,p2;
! 186: UP2 gcd;
! 187:
! 188: p1 = (GF2N)ARG0(arg); p2 = (GF2N)ARG1(arg);
! 189: if ( !p1 )
! 190: *rp = p2;
! 191: else if ( !p2 )
! 192: *rp = p1;
! 193: else {
! 194: gcdup2(p1->body,p2->body,&gcd);
! 195: MKGF2N(gcd,*rp);
1.1 noro 196: }
197: }
198:
199: void Ptdiv(arg,rp)
200: NODE arg;
201: P *rp;
202: {
203: P p1,p2,q1,q2,q,c1,c2,c;
1.5 noro 204: int m;
1.1 noro 205:
206: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg);
207: asir_assert(p1,O_P,"tdiv");
208: asir_assert(p2,O_P,"tdiv");
209: if ( !p1 || !p2 )
210: *rp = 0;
211: else if ( (OID(p1) > O_P) || (OID(p2) > O_P ) )
212: *rp = 0;
1.5 noro 213: else if ( argc(arg) == 3 ) {
214: m = QTOS((Q)ARG2(arg));
215: ptomp(m,p1,&q1); ptomp(m,p2,&q2);
216: if ( divtmp(CO,m,q1,q2,&q) )
217: mptop(q,rp);
218: else
219: *rp = 0;
1.6 noro 220: } else if ( qpcheck((Obj)p1) && qpcheck((Obj)p2) ) {
1.1 noro 221: ptozp(p1,1,(Q *)&c1,&q1); ptozp(p2,1,(Q *)&c2,&q2);
222: if ( divtpz(CO,q1,q2,&q) ) {
223: divq((Q)c1,(Q)c2,(Q *)&c); mulp(CO,q,c,rp);
224: } else
1.6 noro 225: *rp = 0;
226: } else {
227: if ( !divtp(CO,p1,p2,rp) )
1.1 noro 228: *rp = 0;
229: }
230: }
231:
232: void Pudiv(arg,rp)
233: NODE arg;
234: LIST *rp;
235: {
236: P q,r,dnd,dvr;
237: NODE n,tn;
238:
239: asir_assert(ARG0(arg),O_P,"udiv");
240: asir_assert(ARG1(arg),O_P,"udiv");
241: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg);
242: udivpz(dnd,dvr,&q,&r);
243: MKNODE(tn,r,0); MKNODE(n,q,tn); MKLIST(*rp,n);
244: }
245:
246: void Psdivm(arg,rp)
247: NODE arg;
248: Obj *rp;
249: {
250: P q,r,dnd,dnd1,dndm,dvr,dvr1,dvrm,t;
251: V v;
252: VL vl;
253: int m;
254:
255: asir_assert(ARG0(arg),O_P,"sdivm");
256: asir_assert(ARG1(arg),O_P,"sdivm");
257: asir_assert(ARG2(arg),O_N,"sdivm");
258: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg); m = QTOS((Q)ARG2(arg));
259: if ( argc(arg) == 4 ) {
260: v = VR((P)ARG3(arg));
261: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
262: reordvar(CO,v,&vl);
263: ptomp(m,dnd1,&dndm); ptomp(m,dvr1,&dvrm);
264: divsrmp(vl,m,dndm,dvrm,&t,&r); mptop(t,&q);
265: restore_mvar(CO,q,v,(P *)rp);
266: } else {
267: ptomp(m,dnd,&dndm); ptomp(m,dvr,&dvrm);
268: divsrmp(CO,m,dndm,dvrm,&t,&r); mptop(t,(P *)rp);
269: }
270: }
271:
272: void Psremm(arg,rp)
273: NODE arg;
274: Obj *rp;
275: {
276: P q,r,dnd,dnd1,dndm,dvr,dvr1,dvrm,t;
277: V v;
278: VL vl;
279: int m;
280:
281: asir_assert(ARG0(arg),O_P,"sremm");
282: asir_assert(ARG1(arg),O_P,"sremm");
283: asir_assert(ARG2(arg),O_N,"sremm");
284: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg); m = QTOS((Q)ARG2(arg));
285: if ( argc(arg) == 4 ) {
286: v = VR((P)ARG3(arg));
287: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
288: reordvar(CO,v,&vl);
289: ptomp(m,dnd1,&dndm); ptomp(m,dvr1,&dvrm);
290: divsrmp(vl,m,dndm,dvrm,&q,&t); mptop(t,&r);
291: restore_mvar(CO,r,v,(P *)rp);
292: } else {
293: ptomp(m,dnd,&dndm); ptomp(m,dvr,&dvrm);
294: divsrmp(CO,m,dndm,dvrm,&q,&t); mptop(t,(P *)rp);
295: }
296: }
297:
298: void Psqrm(arg,rp)
299: NODE arg;
300: LIST *rp;
301: {
302: P q,q1,r,r1,dnd,dnd1,dndm,dvr,dvr1,dvrm;
303: NODE n,tn;
304: V v;
305: VL vl;
306: int m;
307:
308: asir_assert(ARG0(arg),O_P,"sqrm");
309: asir_assert(ARG1(arg),O_P,"sqrm");
310: asir_assert(ARG2(arg),O_N,"sqrm");
311: dnd = (P)ARG0(arg); dvr = (P)ARG1(arg); m = QTOS((Q)ARG2(arg));
312: if ( argc(arg) == 4 ) {
313: v = VR((P)ARG3(arg));
314: change_mvar(CO,dnd,v,&dnd1); change_mvar(CO,dvr,v,&dvr1);
315: reordvar(CO,v,&vl);
316: ptomp(m,dnd1,&dndm); ptomp(m,dvr1,&dvrm);
1.2 noro 317: divsrmp(vl,m,dndm,dvrm,&q,&r); mptop(q,&q1); mptop(r,&r1);
1.1 noro 318: restore_mvar(CO,q1,v,&q); restore_mvar(CO,r1,v,&r);
319: } else {
320: ptomp(m,dnd,&dndm); ptomp(m,dvr,&dvrm);
1.2 noro 321: divsrmp(CO,m,dndm,dvrm,&q1,&r1); mptop(q1,&q); mptop(r1,&r);
1.1 noro 322: }
323: MKNODE(tn,r,0); MKNODE(n,q,tn); MKLIST(*rp,n);
324: }
325:
326: void Pinva_mod(arg,rp)
327: NODE arg;
328: P *rp;
329: {
330: P dp,f;
331: Q q;
332: int n,i;
333: int mod;
334: V v;
335: UM wf,wdp,winv;
336:
337: asir_assert(ARG0(arg),O_P,"gcda_mod");
338: asir_assert(ARG1(arg),O_N,"gcda_mod");
339: asir_assert(ARG2(arg),O_P,"gcda_mod");
340: dp = (P)ARG0(arg);
341: mod = QTOS((Q)ARG1(arg));
342: f = (P)ARG2(arg);
343: if ( NUM(f) ) {
344: i = invm(rem(NM((Q)f),mod),mod);
345: STOQ(i,q); *rp = (P)q;
346: } else {
347: v = VR(dp);
348: n = MAX(UDEG(dp),UDEG(f));
349: wf = W_UMALLOC(n); wdp = W_UMALLOC(n);
350: winv = W_UMALLOC(n);
351: ptoum(mod,f,wf); ptoum(mod,dp,wdp);
352: invum(mod,wdp,wf,winv);
353: if ( DEG(winv) < 0 )
354: *rp = 0;
355: else {
356: umtop(v,winv,rp);
357: }
358: }
359: }
360:
361: void Psrem_mod(arg,rp)
362: NODE arg;
363: P *rp;
364: {
365: P p1,p2;
366: int n,dr;
367: int mod;
368: V v;
369: UM wp1,wp2,q;
370:
371: asir_assert(ARG0(arg),O_P,"srem_mod");
372: asir_assert(ARG1(arg),O_P,"srem_mod");
373: asir_assert(ARG2(arg),O_N,"srem_mod");
374: p1 = (P)ARG0(arg); p2 = (P)ARG1(arg); mod = QTOS((Q)ARG2(arg));
375: if ( !p1 || NUM(p1) )
376: *rp = p1;
377: else {
378: v = VR(p1);
379: n = MAX(UDEG(p1),UDEG(p2));
380: wp1 = W_UMALLOC(n); wp2 = W_UMALLOC(n); q = W_UMALLOC(n);
381: ptoum(mod,p1,wp1); ptoum(mod,p2,wp2);
382: dr = divum(mod,wp1,wp2,q);
383: if ( ( DEG(wp1) = dr ) == -1 )
384: *rp = 0;
385: else
386: umtop(v,wp1,rp);
387: }
388: }
389:
390: void Purem(arg,rp)
391: NODE arg;
392: P *rp;
393: {
394: asir_assert(ARG0(arg),O_P,"urem");
395: asir_assert(ARG1(arg),O_P,"urem");
396: uremp((P)ARG0(arg),(P)ARG1(arg),rp);
397: }
398:
399: void Pugcd(arg,rp)
400: NODE arg;
401: P *rp;
402: {
403: asir_assert(ARG0(arg),O_P,"ugcd");
404: asir_assert(ARG1(arg),O_P,"ugcd");
405: ugcdp((P)ARG0(arg),(P)ARG1(arg),rp);
406: }
407:
408: void invum(mod,dp,f,inv)
409: int mod;
410: UM dp,f,inv;
411: {
412: UM g1,g2,a1,a2,a3,wm,q,tum;
413: int d,dr;
414:
415: d = DEG(dp)+DEG(f)+10;
416: g1 = W_UMALLOC(d); g2 = W_UMALLOC(d); a1 = W_UMALLOC(d);
417: a2 = W_UMALLOC(d); a3 = W_UMALLOC(d); wm = W_UMALLOC(d);
418: q = W_UMALLOC(d);
419: DEG(a1) = 0; COEF(a1)[0] = 1; DEG(a2) = -1;
420: cpyum(f,g1); cpyum(dp,g2);
421: while ( 1 ) {
422: dr = divum(mod,g1,g2,q); tum = g1; g1 = g2; g2 = tum;
423: if ( ( DEG(g2) = dr ) == -1 )
424: break;
425: mulum(mod,a2,q,wm); subum(mod,a1,wm,a3); dr = divum(mod,a3,dp,q);
426: tum = a1; a1 = a2; a2 = a3; a3 = tum; DEG(a3) = dr;
427: }
428: if ( DEG(g1) != 0 )
429: DEG(inv) = -1;
430: else if ( COEF(g1)[0] != 1 )
431: mulsum(mod,a2,invm(COEF(g1)[0],mod),inv);
432: else
433: cpyum(a2,inv);
434: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>
![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to pdiv.c CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib2 / asir2000 / builtin